Juli 2015

Band 30, Nummer 7

Programmiersprache R – Einführung in R für C#-Programmierer

Von James McCaffrey

Die Sprache R wird von Data Scientists und Programmierern für statistische Berechnungen verwendet. Teilweise aufgrund der von Softwaresystemen erfassten wachsenden Datenmengen und der Notwendigkeit, diese Daten zu analysieren, ist R eine der am schnellsten wachsenden Technologien unter meinen Kollegen, die C# verwenden. Vertrautheit mit R kann eine wertvolle Ergänzung Ihrer technischen Qualifikationen sein.

Die Sprache R ist ein GNU-Projekt und stellt kostenlose Software dar. R wurde von einer Sprache mit dem Namen S (für „Statistik“) abgeleitet, die in den Bell Laboratories in den 1970er Jahren entwickelt wurde. Es gibt viele hervorragende Onlinetutorials für R, aber die meisten dieser Tutorials gehen davon aus, dass Sie an einer Universität Statistik studieren. Dieser Artikel hat das Ziel, C#-Programmierern zu helfen, so schnell wie möglich einen flüssigen Umgang mit R zu erwerben.

Die Zielrichtung dieses Artikels können Sie am leichtesten erfassen, indem Sie sich die R-Beispielsitzung in Abbildung 1 ansehen. Die Sitzung hat zwei nicht miteinander zusammenhängende Themen. Die ersten paar Befehlssätze zeigen einen sogenannten Chi-Quadrat-Test (im Text als „chi-square“ oder auch als „chi-squared“ bezeichnet) für einheitliche Verteilung. Der zweite Satz von Befehlen zeigt ein Beispiel für lineare Regression, die in meinen Augen im Computing das „Hello World“ der Statistik darstellt.

R-Beispielsitzung
Abbildung 1 R-Beispielsitzung

Die R-Website finden Sie unter r-project.org. Die Website weist Links zu verschiedenen Spiegelwebsites auf, von denen Sie R herunterladen und installieren können. Das Installationspaket besteht aus einer einfachen, selbstextrahierenden Programmdatei. R wird offiziell unter Windows XP und höher unterstützt, ebenso wie auf den meisten anderen gebräuchlichen (nicht Windows-) Plattformen. Ich habe R auf Computern mit Windows 7 und Windows 8 ohne Probleme installiert. Standardmäßig stellt der Installationsprozess sowohl die 32-Bit- als auch die 64-Bit-Version bereit.

In diesem Artikel wird davon ausgegangen, dass Sie zumindest fortgeschrittene C#-Programmierkenntnisse haben (sodass Sie die Erläuterungen der Ähnlichkeiten und Unterschiede zwischen C# und R nachvollziehen können), aber es wird nicht angenommen, dass Sie irgendetwas über R wissen. Ein C#-Demoprogramm wird zur Gänze vorgestellt, das auch in der Downloaddatei, die diesen Artikel begleitet, verfügbar ist.

Der Chi-Quadrat-Test in R

Der erste Blick auf Abbildung 1 offenbart, dass die Verwendung von R sich deutlich von der Verwendung von C# unterscheidet. Es ist zwar möglich, Skripts in R zu erstellen, jedoch wird R meistens im interaktiven Modus in einer Befehlsshell verwendet. Das erste R-Beispiel ist eine Analyse zur Prüfung, ob ein normaler sechsseitiger Würfel in Ordnung ist oder nicht. Bei vielen Würfen sollte ein ordnungsgemäß hergestellter Würfel der Annahme entsprechend ungefähr die gleiche Instanzanzahl für jedes der sechs möglichen Ergebnisse erreichen.

Die Eingabeaufforderung von R wird durch das Symbol (>) in der Shell angezeigt. Die erste in Abbildung 1 eingegebene Anweisung beginnt mit dem Zeichen (#), das in R als Symbol zum Kennzeichnen eines Kommentars dient.

Der erste eigentliche Befehl in der Beispielsitzung ist:

> observed <- c(20, 28, 12, 32, 22, 36)

Dadurch wird mithilfe der Funktion „c“ (für concatenate) ein Vektor mit dem Namen „observed“ erstellt. Ein Vektor enthält Objekte vom gleichen Datentyp. Eine ungefähr gleichbedeutende C#-Anweisung wäre:

var observed = new int[] { 20, 28, 12, 32, 22, 36 };

Der erste Wert, 20, ist die Anzahl der Vorkommen der Punktzahl 1, 28 ist die Anzahl der Vorkommen der Punktzahl 2 usw. Die Summe der sechs Anzahlen von Vorkommen ist 20 + 28 + . . + 36 = 150. Man sollte erwarten, dass ein ordnungsgemäß gefertigter Würfel eine Anzahl von ungefähr 150/6 = 25 für jedes mögliche Ergebnis ergibt. Aber die Anzahl der beobachteten Würfe mit Punktzahl 3 (12) sieht verdächtig niedrig aus.

Nach dem Erstellen des Vektors wird ein Chi-Quadrat-Test mithilfe der Funktion „chisq.test“ durchgeführt:

> chisq.test(observed)

In R wird häufig der Punkt (.) anstelle des Unterstrichzeichens (_) verwendet, um einfacher lesbare Variablen- und Funktionsnamen zu erstellen. Das Ergebnis des Aufrufs der Funktion „chisq.test“ ist eine ganze Menge Text:

Chi-Quadrat-Test für die bestehenden Wahrscheinlichkeiten

data:  observed
X-squared = 15.28, df = 5, p-value = 0.009231

In der Terminologie von C# ausgedrückt geben die meisten R-Funktionen eine Datenstruktur zurück, die ignoriert werden kann, enthalten aber auch viele Anweisungen, die zu „Console.WriteLine“ äquivalent sind und Ausgaben verfügbar machen. Beachten Sie, dass es Ihre Aufgabe ist, die Bedeutung von R-Ausgaben zu entschlüsseln. Weiter unten in diesem Artikel demonstriere ich die Erstellung des gleichbedeutenden Chi-Quadrat-Tests in reinem C#-Code (ohne Bibliotheken).

In diesem Beispiel ist der „X-squared“-Wert von 15,28 die berechnete Chi-Quadrat-Statistik (der griechische Buchstabe „chi“ ist dem großen „X“ ähnlich). Der Wert „0,0“ zeigt an, dass die beobachteten Werte genau dem erwarteten Ergebnis bei einem einwandfreien Würfel entsprechen. Größere Werte von Chi-Quadrat zeigen eine zunehmende Wahrscheinlichkeit an, dass die beobachteten Fallzahlen bei einem einwandfreien Würfel nicht zufällig zustande gekommen sind. Der df-Wert 5 bezeichnet die Freiheitgrade („degrees of freedom“), einen weniger als die Anzahl der beobachteten Werte. Der df-Wert ist für diese Analyse nicht allzu wichtig.

Der p-Wert von 0,009231 ist die Wahrscheinlichkeit, dass die beobachteten Fallzahlen bei einem einwandfreien Würfel zufällig zustande gekommen sind. Da der p-Wert so klein ist, unter 1 Prozent, würde man schließen, dass die beobachteten Werte höchst wahrscheinlich nicht zufällig zustande gekommen sind und es daher einen statistischen Beleg dafür gibt, dass der fragliche Würfel eine Neigung hat.

Lineare Regressionsanalyse mit R

Der zweite Satz Anweisungen in Abbildung 1 zeigt ein Beispiel für lineare Regression. Lineare Regression ist eine statistische Technik, die zum Beschreiben der Beziehung zwischen einer numerischen Variablen (in der Statistik als „abhängige Variable“ bezeichnet) und einer oder mehreren erklärenden Variablen (die man unabhängige Variablen nennt), die entweder numerisch oder kategorisch sein können, verwendet wird. Wenn nur eine unabhängige erklärende/Prädiktorvariable vorhanden ist, wird die Technik als einfache lineare Regression bezeichnet. Bei zwei oder mehr unabhängigen Variablen, wie im Beispiel der Demo, wird die Technik als multiple lineare Regression bezeichnet.

Vor dem Ausführen der linearen Regressionsanalyse habe ich eine aus acht Elementen bestehende, kommagetrennte Textdatei mit dem Namen „DummyData.txt“ im Verzeichnis „C:\IntroToR“ erstellt, die diesen Inhalt aufweist:

Color,Length,Width,Rate
blue, 5.4, 1.8, 0.9
blue, 4.8, 1.5, 0.7
blue, 4.9, 1.6, 0.8
pink, 5.0, 1.9, 0.4
pink, 5.2, 1.5, 0.3
pink, 4.7, 1.9, 0.4
teal, 3.7, 2.2, 1.4
teal, 4.2, 1.9, 1.2

Die Datei soll Blumendaten mit der Farbe der Blume, der Länge und Breite des Blütenblatts und der Wachstumsrate darstellen. Dem liegt die Absicht zugrunde, den Wert für die Wachstumsrate (in der letzten Spalte) aus den Werten für Farbe, Länge und Breite vorherzusagen. Nach einer Kommentaranweisung sind die ersten drei R-Befehle in der linearen Regressionsanalyse:

> setwd("C:\\IntroToR")
> data <- read.table("DummyData.txt",
  header=TRUE, sep=",")
> print(data)

Der erste Befehl legt das Arbeitsverzeichnis fest, damit der Pfad zur Quelldatendatei nicht voll qualifiziert werden muss. Anstelle des Symbols (\\), das in C# gebräuchlich ist, hätte ich (/) verwenden können, das auf anderen Plattformen als Windows gängig ist.

Der zweite Befehl lädt die Daten in den Arbeitsspeicher in ein Tabellenobjekt mit dem Namen „data“. Beachten Sie, dass R benannte Parameter verwendet. Der Headerparameter gibt an, ob die erste Zeile aus Headerinformationen besteht („TRUE“ oder kurz „T“) oder nicht („FALSE“ oder „F“). R unterscheidet Groß- und Kleinschreibung, und Sie können normalerweise entweder den Operator (<-) oder den Operator (=) zum Zuweisen von Werten verwenden. Das bleibt im Wesentlichen Ihrem Geschmack überlassen. Ich verwende meistens (<-) für die Zuweisung von Objekten und (=) für die Zuweisung von Parameterwerten.

Der Parameter „sep“ (separator) gibt an, wie die Werte in den einzelnen Zeilen getrennt sind. Beispielsweise würde (\t) Werte mit Tabulatortrennung und (" ") Werte mit Trennung durch Leerzeichen angeben.

Die Funktion „print“ zeigt die Datentabelle im Arbeitsspeicher an. Die Funktion „print“ weist viele optionale Parameter auf. Beachten Sie, dass die Ausgabe in Abbildung 1 Datenelementindizes beginnend mit 1 anzeigt. Für Array-, Matrix- und Objektindizes ist R eine 1-basierte Sprache, nicht etwa 0-basiert, wie C#.

Die lineare Regressionsanalyse wird mit diesen zwei R-Befehlen durchgeführt:

> model <- lm(data$Rate ~ (data$Color + data$Length + data$Width))
> summary(model)

Der erste Befehl lässt sich so lesen: „Speichere in einem Objekt mit dem Namen ‚model‘ das Ergebnis der lm-Funktionsanalyse (lineares Modell), bei der die vorherzusagende abhängige Variable die Spalte ‚Rate‘ im Tabellenobjekt (data$Rate) darstellt und die unabhängigen Prädiktorvariablen ‚Color‘, ‚Length‘ und ‚Width‘ sind“. Der zweite Befehl bedeutet „Zeige nur die einfachen Ergebnisse der im Objekt mit dem Namen ‚model‘ gespeicherten Analyse an“.

Die lm-Funktion generiert große Mengen an Informationen. Angenommen, Sie möchten den Wert für „Rate“ vorhersagen, wenn die Eingabewerte „Color = pink“, „Length = 5.0“ und „Width = 1.9“ sind. (Beachten Sie, dass dies dem Datenelement [4] entspricht, das einen tatsächlichen Wert für „Rate“ von 0,4 aufweist.) Zum Treffen einer Vorhersage verwenden Sie dann die Werte in der Spalte „Estimate“:

Coefficients:
               Estimate Std. Error t value Pr(>|t|)
(Intercept)    -0.14758    0.48286  -0.306  0.77986
data$Colorpink -0.49083    0.04507 -10.891  0.00166 **
data$Colorteal  0.35672    0.09990   3.571  0.03754 *
data$Length     0.04159    0.07876   0.528  0.63406
data$Width      0.45200    0.11973   3.775  0.03255 *

Wenn X für die unabhängigen Variablenwerte und Y für die vorhergesagte Rate steht, dann gilt:

X = (blue = NA, pink = 1, teal = 0, Length = 5.0, Width = 1.9)
Y = -0.14758 + (-0.49083)(1) + (0.35672)(0) + (0.04159)(5.0) + (0.45200)(1.9)
  = -0.14758 + (-0.49083) + (0) + (0.20795) + (0.85880)
  = 0.42834

Beachten Sie, dass die vorhergesagte Rate 0,43 recht nahe an der tatsächlichen Rate 0,40 liegt.

In Worten ausgedrückt berechnen Sie zum Treffen einer Vorhersage mithilfe des Modells eine lineare Summe von Produkten der Werte von „Estimate“, die mit ihren entsprechenden X-Werten multipliziert wurden. Der Wert von „Intercept“ ist eine Konstante, die keiner Variablen zugeordnet ist. Wenn Sie kategorische erklärende Variablen verwenden, wird einer der Werte fortgelassen (in diesem Fall Blau).

Die Informationen unten im Ausgabedisplay zeigen an, wie gut die unabhängigen Variablen „Color“, „Length“ und „Width“ die abhängige Variable erklären, Rate:

Residual standard error: 0.05179 on 3 degrees of freedom
Multiple R-squared:  0.9927,    Adjusted R-squared:  0.9829
F-statistic: 101.6 on 4 and 3 DF,  p-value: 0.00156

Der Wert „Multiple R-squared“ (0,9927) ist der Prozentsatz der Varianz in der abhängigen Variablen, die sich durch die lineare Kombination der unabhängigen Variablen erklärt. Etwas anders ausgedrückt stellt „R-squared“ einen Wert zwischen 0 und 1 dar, bei dem höhere Werte ein besseres Vorhersagemodell bedeuten. In diesem Fall ist der Wert für „R-squared“ extrem hoch, was anzeigt, dass „Color“, „Length“ und „Width“ den Wert für „Rate“ sehr genau vorhersagen können. Die Werte „F-statistic“, „Adjusted R-squared“ und „p-value“ sind weitere Kennzahlen der Modellgüte.

Eine Absicht dieses Beispiels ist, zu zeigen, dass Ihre mit Abstand größte Herausforderung beim Programmieren mit R im Verstehen der Statistik hinter den Funktionen der Sprache liegt. Die meisten Menschen lernen R inkrementell, indem sie jeweils das Wissen zu einer Technik neu erwerben, wie es zum Beantworten einer bestimmten Frage erforderlich ist. Ein analoges Verhalten in C# wäre, die verschiedenen Sammlungsobjekte im Namespace „Collections.Generic“ zu lernen, etwa die Hashtable- und Queue-Klassen. Die meisten Entwickler lernen immer jeweils eine Datenstruktur, statt Informationen über alle Klassen zu memorieren, bevor sie eine von ihnen in einem Projekt verwenden.

Ein weiterer Chi-Quadrat-Test

Die Art von Chi-Quadrat-Test in Abbildung 1 wird oft als Test auf gleichmäßige Verteilung bezeichnet, weil er überprüft, ob die beobachteten Daten alle gleiche Häufigkeiten aufweisen; also, ob die Daten gleichmäßig verteilt sind. Es gibt eine Reihe weiterer Chi-Quadrat-Tests, darunter einer, der als Chi-Quadrat-Unabhängigkeitstest bezeichnet wird.

Angenommen, Sie untersuchen eine Gruppe von 100 Personen und möchten wissen, ob das Geschlecht (männlich, weiblich) von der Nähe zu politischen Parteien (Demokraten, Republikaner, Andere) unabhängig ist. Stellen Sie sich Ihre Daten als Kontingenztafel vor, wie in Abbildung 2. Offenbar scheinen Männer eher den Republikanern zugeneigt zu sein als Frauen.

Abbildung 2 Kontingenztafel

  Dem Rep Andere  
Männlich 15 25 10 50
Weiblich 30 15 5 50
  45 40 15 100

Um R für den Test zu verwenden, ob die zwei Faktoren, Geschlecht und Zugehörigkeit, statistisch unabhängig sind, platzieren Sie im ersten Schritt die Daten mit diesem Befehl in einer numerischen Matrix:

> cm <- matrix( c(15,30,25,15,10,5), nrow=2, ncol=3 )

Beachten Sie, dass Matrixdaten in R nach Spalten (von oben nach unten, von links nach rechts) statt nach Zeilen (von links nach rechts, von oben nach unten) gespeichert sind, wie in C#.

Der R-Befehl für den Chi-Quadrat-Test lautet:

> chisq.test(cm)

Das Ergebnis für den p-Wert des Tests ist 0,01022, was bei einem Signifikanzniveau von 5 Prozent darauf hinweist, dass die beiden Faktoren nicht unabhängig voneinander sind. Anders gesagt, es besteht eine statistische Beziehung zwischen Geschlecht und Parteizugehörigkeit.

Beachten Sie, dass im ersten Chi-Quadrat-Beispiel mit dem Würfel der Eingabeparameter ein Vektor ist, im zweiten Beispiel mit Geschlecht und Zugehörigkeit die Eingabe aber aus einer Matrize besteht. Die Funktion „chisq.test“ weist insgesamt sieben Parameter auf, von denen einer erforderlich ist (ein Vektor oder eine Matrize), gefolgt von sechs optionalen benannten Parametern.

Wie viele C#-Methoden in den Microsoft .NET Framework-Namespaces, sind die meisten R-Funktionen stark überladen. In C# wird Überladung normalerweise mithilfe mehrerer Methoden mit gleichen Namen, aber verschiedenen Parametern implementiert. Die Verwendung von Generika ist auch eine Form der Überladung. In R wird Überladung in Form einer einzelnen Funktion mit vielen optionalen, benannten Parametern implementiert.

Diagramme mit R

Wenn ich als C#-Programmierer ein Diagramm aus Ausgabedaten erstellen möchte, führe ich normalerweise mein Programm aus, kopiere die Ausgabedaten, füge sie mit STRG+V in den Editor ein, um überflüssige Steuerzeichen zu entfernen, kopiere diese Daten, füge sie in Excel ein, und erstelle dann ein Diagramm mithilfe von Excel. Das ist ein bisschen gehackt, aber es funktioniert in den meisten Situationen gut.

Eine der Stärken des R-Systems ist eine systemeigene Fähigkeit zum Erstellen von Diagrammen. Sehen Sie sich das Beispiel in Abbildung 3 an. Diese Art Diagramm ist in Excel nicht ohne ein geeignetes Add-In möglich.

3D-Diagramm mit R
Abbildung 3 3D-Diagramm mit R

Über das in Abbildung 1 gezeigte Shellprogramm hinaus verfügt R außerdem über eine Halb-GUI-Benutzeroberfläche, „RGui.exe“, die zum Erstellen von Diagrammen verwendet wird.

Das Diagramm in Abbildung 3 bildet die Funktion „z = f(x,y) = x * e^(-(x^2 + y^2))“ ab. Die ersten drei R-Befehle zum Erstellen des Diagramms sind:

> rm(list=ls())
> x <- seq(-2, 2, length=25)
> y <- seq(-2, 2, length=25)

Die rm-Funktion löscht ein Objekt aus dem aktuellen Arbeitsbereich im Arbeitsspeicher. Der verwendete Befehl ist eine magische R-Beschwörung zum Löschen aller Objekte. Der zweite und der dritte Befehl erstellen Vektoren aus 25 Werten mit gleichen Abständen von -2 bis +2 einschließlich. Stattdessen hätte auch die c-Funktion verwendet werden können.

Die nächsten beiden Befehle sind:

> f <- function(x,y) { x * exp(-(x^2
  + y^2)) }
> z <- outer(x,y,f)

Der erste Befehl zeigt die Definition einer Funktion in R mithilfe des Schlüsselworts „function“. Die integrierte R-Funktion mit dem Namen „outer“ erstellt eine Wertmatrix mithilfe der Vektoren „x“ und „y“ und einer Funktionsdefinition „f“. Das Ergebnis ist eine 25 x 25-Matrix, bei der der Wert in jeder Zelle der Wert der Funktion „f“ ist, der „x“ und „y“ entspricht.

Die nächsten beiden Befehle sind:

> nrz <- nrow(z)
> ncz <- ncol(z)

Die Funktionen „nrow“ und „ncol“ geben die Anzahl der Zeilen oder die Anzahl der Spalten in ihrem Matrixargument zurück. In diesem Fall sind beide Werte 25.

Der nächste R-Befehl verwendet die Funktion „colorRampPalette“, um eine benutzerdefinierte Farbverlaufspalette zum Zeichnen des Diagramms zu erstellen:

> jet.colors <- colorRampPalette(c("midnightblue", "blue",
+ "cyan", "green", "yellow", "orange", "red", "darkred"))

Wenn Sie eine lange Zeile in R eingeben, springt der Cursor des Systems beim Eingeben der <EINGABETASTE> nach unten in die nächste Zeile und platziert ein +-Zeichen als Eingabeaufforderung, um anzuzeigen, dass der Befehl noch nicht abgeschlossen ist. Weiter:

> nbcol <- 64
> color <- jet.colors(nbcol)

Das Ergebnis dieser beiden Befehle ist ein Vektor mit dem Namen „color“, der 64 verschiedene Farbwerte enthält, die von einem sehr dunklen Blau über Grün und Gelb bis zu einem sehr dunklen Rot reichen. Weiter:

> zfacet <- z[-1,-1] + z[-1,-ncz] + z[-nrz,-1] + z[-nrz,-ncz]
> facetcol <- cut(zfacet,nbcol)

Wenn Sie das Diagramm in Abbildung 3 genauer betrachten, sehen Sie, dass die Fläche aus 25 x 25 = 625 kleinen Quadraten oder Facetten besteht. Die beiden vorhergehenden Befehle erstellen eine 25 x 25-Matrix, bei der der Wert jeder Zelle eine von 64 Farben ist, die für die entsprechende Facette der Oberfläche verwendet wird.

Schließlich wird das 3D-Diagramm mithilfe der Funktion „persp“ (perspektivisches Diagramm) angezeigt:

> persp(x,y,z, col=color[facetcol], phi=20, theta=-35,
+ ticktype="detailed", d=5, r=1, shade=0.1, expand=0.7)

Die Funktion „persp“ weist viele optionale, benannte Parameter auf. Der Parameter „col“ gibt die zu verwendende(n) Farbe oder Farben an. Die Parameter „phi“ und „theta“ legen den Betrachtungswinkel (links und rechts sowie oben und unten) des Diagramms fest. Der Parameter „ticktype“ steuert, wie die Werte auf den x-, y- und z-Achsen angezeigt werden. Die Parameter „r“ und „d“ steuern die wahrgenommene Augendistanz zum Diagramm und den wahrgenommenen 3D-Effekt. Der Parameter mit dem Namen „shade“ steuert die simulierte Abschattung von einer virtuellen Lichtquelle. Der Parameter „expand“ steuert das Verhältnis von Höhe und Breite des Diagramms. Die Funktion „persp“ weist noch viele weitere Parameter auf, aber die hier verwendeten reichen für die meisten Situationen aus.

Dieses Beispiel macht deutlich, dass R über sehr leistungsfähige und flexible systemeigene Möglichkeiten zur Diagrammerstellung verfügt, sie befinden sich jedoch auf ziemlich niederer Ebene und erfordern erheblichen Aufwand.

Der Chi-Quadrat-Test in C#

Um Verständnis für die Ähnlichkeiten und Unterschiede zwischen Analysen in der Sprache R und Programmierung in C# zu entwickeln, ist es nützlich, eine C#-Implementierung eines Chi-Quadrat-Tests zu betrachten. Darüber hinaus kann der C#-Code eine willkommene Erweiterung Ihrer persönlichen Codebibliothek darstellen. Sehen Sie sich das C#-Demoprogramm in Abbildung 4 an.

Chi-Quadrat-Test mithilfe von C#
Abbildung 4 Chi-Quadrat-Test mithilfe von C#

Das Demoprogramm stellt eine Annäherung an den Chi-Quadrat-Test in der Sprache R dar, der in Abbildung 1 dargestellt ist. Beachten Sie, dass die Ausgabewerte der C#-Demo exakt die gleichen wie in der R-Sitzung sind.

Zum Erstellen des Demoprogramms wurde in Visual Studio ein neues C#-Konsolenanwendungprojekt mit dem Namen „ChiSquare“ erstellt. Nach dem Laden des Vorlagencodes in den Editor habe ich im Fenster des Projektmappen-Explorers mit der rechten Maustaste auf „Program.cs“ geklickt, die Datei in „ChiSquareProgram.cs“ umbenannt und Visual Studio die automatische Umbenennung der Klasse „Program“ in „ChiSquareProgram“ gestattet.

Das vollständige C#-Demoprogramm ist in Abbildung 5 aufgelistet. Sie werden feststellen, dass der Quellcode länger ist, als man erwartet. Das Implementieren von statistischer Programmierung in reinem C# ist nicht überwältigend schwierig, der Code neigt aber zur Länge.

Abbildung 5 C# Chi-Quadrat-Demoprogramm

using System;
namespace ChiSquare
{
  class ChiSquareProgram
  {
    static void Main(string[] args)
    {
      try
      {
        Console.WriteLine("\nBegin Chi-square test using C# demo\n");
        Console.WriteLine(
          "Goal is to see if one die from a set of dice is biased or not\n");
        int[] observed = new int[] { 20, 28, 12, 32, 22, 36 };
        Console.WriteLine("\nStarting chi-square test");
        double p = ChiSquareTest(observed);
        Console.WriteLine("\nChi-square test complete");
        double crit = 0.05;
        if (p < crit)
        {
          Console.WriteLine("\nBecause p-value is below critical value of " +
            crit.ToString("F2"));
          Console.WriteLine("the null hypothsis is rejected and we conclude");
          Console.WriteLine("the data is unlikely to have happened by chance.");
        }
        else
        {
          Console.WriteLine("\nBecause p-value is not below critical value of " +
            crit.ToString("F2"));
          Console.WriteLine(
            "the null hypothsis is accepted (not rejected) and we conclude");
          Console.WriteLine("the observed data could have happened by chance.");
        }
        Console.WriteLine("\nEnd\n");
        Console.ReadLine();
      }
      catch (Exception ex)
      {
        Console.WriteLine(ex.Message);
        Console.ReadLine();
      }
    } // Main
    static void ShowVector(int[] v)
    {
      for (int i = 0; i < v.Length; ++i)
        Console.Write(v[i] + " ");
      Console.WriteLine("");
    }
    static double ChiSquareTest(int[] observed)
    {
      Console.WriteLine("Observed frequencies are: ");
      ShowVector(observed);
      double x = ChiSquareStatistic(observed);
      Console.WriteLine("\nX-squared = " + x.ToString("F2"));
      int df = observed.Length - 1;
      Console.WriteLine("\ndf = " + df);
      double p = ChiSquareProb(x, df);
      Console.WriteLine("\np-value = " + p.ToString("F6"));
      return p;
    }
    static double ChiSquareStatistic(int[] observed)
    {
      double sumObs = 0.0;
      for (int i = 0; i < observed.Length; ++i)
        sumObs += observed[i];
      double expected = (int)(sumObs / observed.Length);
      double result = 0.0;
      for (int i = 0; i < observed.Length; ++i)
      {
        result += ((observed[i] - expected) *
         (observed[i] - expected)) / expected;
      }
      return result;
    }
    public static double ChiSquareProb(double x, int df)
    {
      // x = a computed chi-square value. df = degrees of freedom.
      // output = prob. the x value occurred by chance.
      // So, for example, if result < 0.05 there is only a 5% chance
      // that the x value occurred by chance and, therefore,
      // we conclude that the actual data which produced x is
      // NOT the same as the expected data.
      // This function can be used to create a ChiSquareTest procedure.
      // ACM Algorithm 299 and update ACM TOMS June 1985.
      // Uses custom Exp() function below.
      if (x <= 0.0 || df < 1)
        throw new Exception("parameter x must be positive " +
        "and parameter df must be 1 or greater in ChiSquaredProb()");
      double a = 0.0; // 299 variable names
      double y = 0.0;
      double s = 0.0;
      double z = 0.0;
      double e = 0.0;
      double c;
      bool even; // Is df even?
      a = 0.5 * x;
      if (df % 2 == 0) even = true; else even = false;
      if (df > 1) y = Exp(-a); // ACM update remark (4)
      if (even == true) s = y; else s = 2.0 * Gauss(-Math.Sqrt(x));
      if (df > 2)
      {
        x = 0.5 * (df - 1.0);
        if (even == true) z = 1.0; else z = 0.5;
        if (a > 40.0) // ACM remark (5)
        {
          if (even == true) e = 0.0;
          else e = 0.5723649429247000870717135;
          c = Math.Log(a); // log base e
          while (z <= x)
          {
            e = Math.Log(z) + e;
            s = s + Exp(c * z - a - e); // ACM update remark (6)
            z = z + 1.0;
          }
          return s;
        } // a > 40.0
        else
        {
          if (even == true) e = 1.0;
          else e = 0.5641895835477562869480795 / Math.Sqrt(a);
          c = 0.0;
          while (z <= x)
          {
            e = e * (a / z); // ACM update remark (7)
            c = c + e;
            z = z + 1.0;
          }
          return c * y + s;
        }
      } // df > 2
      else
      {
        return s;
      }
    } // ChiSquare()
    private static double Exp(double x) // ACM update remark (3)
    {
      if (x < -40.0) // ACM update remark (8)
        return 0.0;
      else
        return Math.Exp(x);
    }
    public static double Gauss(double z)
    {
      // input = z-value (-inf to +inf)
      // output = p under Normal curve from -inf to z
      // e.g., if z = 0.0, function returns 0.5000
      // ACM Algorithm #209
      double y; // 209 scratch variable
      double p; // result. called 'z' in 209
      double w; // 209 scratch variable
      if (z == 0.0)
        p = 0.0;
      else
      {
        y = Math.Abs(z) / 2;
        if (y >= 3.0)
        {
          p = 1.0;
        }
        else if (y < 1.0)
        {
          w = y * y;
          p = ((((((((0.000124818987 * w
            - 0.001075204047) * w + 0.005198775019) * w
            - 0.019198292004) * w + 0.059054035642) * w
            - 0.151968751364) * w + 0.319152932694) * w
            - 0.531923007300) * w + 0.797884560593) * y * 2.0;
        }
        else
        {
          y = y - 2.0;
          p = (((((((((((((-0.000045255659 * y
            + 0.000152529290) * y - 0.000019538132) * y
            - 0.000676904986) * y + 0.001390604284) * y
            - 0.000794620820) * y - 0.002034254874) * y
            + 0.006549791214) * y - 0.010557625006) * y
            + 0.011630447319) * y - 0.009279453341) * y
            + 0.005353579108) * y - 0.002141268741) * y
            + 0.000535310849) * y + 0.999936657524;
        }
      }
      if (z > 0.0)
        return (p + 1.0) / 2;
      else
        return (1.0 - p) / 2;
    } // Gauss()
  } // Program class
} // ns

Die Methode „Main“ besteht hauptsächlich aus „WriteLine“-Anweisungen. Der eigentliche aufrufende Code ist:

int[] observed = new int[] { 20, 28, 12, 32, 22, 36 };
double p = ChiSquareTest(observed);
double crit = 0.05;
if (p < crit) {
  // Messages
} else {
  // Messages
}

Der C#-Chi-Quadrat-Test akzeptiert ein Array von beobachteten Werten, berechnet den statistischen Chi-Quadrat-Wert und zeigt ihn an; er berechnet die Freiheitsgrade und zeigt sie an; und er berechnet den p-Wert und gibt ihn zurück.

Die Methode „ChiSquareTest“ ruft drei Hilfsmethoden auf:

ShowVector(observed);
double x = ChiSquareStatistic(observed);
int df = observed.Length - 1;
double p = ChiSquareProb(x, df);

Die Methode „ShowVector“ zeigt den Eingabevektor an, ähnlich dem in der Funktion „chisq.test“ in R verwendeten Ansatz zur Ausgabe (echoing) der Werte der Eingabeparameter. Die Methode „ChiSquareStatistic“ gibt das berechnete Chi-Quadrat zurück („X-squared“ in der R-Ausgabe), und die Methode „ChiSquare­Prob“ verwendet den Rückgabewert von „ChiSquareStatistic“, um eine Wahrscheinlichkeit zu berechnen (den „p-value“ der R-Ausgabe).

Die Methode „ChiSquareStatistic“ ist ein einfacher Test auf gleichmäßige Verteilung. Die Statistikformel für die Chi-Quadrat-Statistik ist:

chi-squared = Sum( (observed - expected)^2 / expected )

Vor dem Berechnen von „chi-squared“ müssen Sie also die erwarteten Werte berechnen, die den beobachteten Werten zugeordnet werden. Wie zuvor erläutert, addieren Sie den Wert „count“ im Array „observed“ (150 in der Demo), dividieren dann die Summe durch die Anzahl der Werte im Array (6 in der Demo), um den erwarteten Wert (25) für alle Zellen zu erhalten.

Der schwierigste Teil beim Erstellen eines Chi-Quadrat-Tests besteht in der Berechnung von „p-value“. Wenn Sie die Definition der Methode „ChiSquare­Prob“ in Abbildung 5 durchsehen, werden Sie schnell feststellen, dass sie tiefgehende Spezialkenntnisse erfordert. Darüber hinaus ruft die Methode eine Hilfsmethode mit den Namen „Gauss“ auf, die gleichermaßen komplex ist.

Das Codieren komplizierter numerischer Funktionen ist tatsächlich ziemlich einfach, weil die meisten Funktionen seit Jahrzehnten gelöst sind. Zwei der von mir besonders häufig verwendeten Ressourcen sind das ACM (Association for Computing Machinery)-Repository der gesammelten Algorithmen und das „Handbook of Mathematical Functions“ von Abramowitz und Stegun (Dover Publications, 1965) – es ist so bekannt, dass es oft nur „A&S“ genannt wird. Beide Referenzen sind kostenlos an mehreren Orten im Web erhältlich.

Beispielsweise implementiert die Methode „ChiSquareProb“ den ACM-Algorithmus Nr. 299, und die Gauss-Hilfsmethode implementiert den ACM-Algorithmus Nr. 209. Eine Alternative zu ACM Nr. 209 ist die A&S-Formel Nr. 7.1.26 zuzüglich einer einfachen Wrappinganweisung.

Viele der Funktionen in den ACM-Algorithmen und A&S sind in bestehenden C#-Bibliotheken implementiert; die meisten dieser Bibliotheken sind allerdings ziemlich umfangreich. Nach Möglichkeit erstelle ich Code von Grund auf und verwende nur die benötigten Methoden, um externe Abhängigkeiten zu vermeiden.

Ein paar Kommentare

Dieser Artikel zeigt nur einen kleinen Bruchteil der Sprache R, aber er sollte ausreichen, um Ihnen den Einstieg zu ermöglichen. C#-Programmierer kommen mit R normalerweise auf eine von zwei Weisen in Berührung. Erstens können Sie R verwenden, um eigene statistische Analysen durchzuführen. Wie dieser Artikel zeigt, ist die Verwendung von R recht einfach, vorausgesetzt, Sie verstehen die zugrundeliegende Statistik. Zweitens kommt es vor, dass Sie mit jemand zusammenarbeiten müssen, der R als primäre Sprache verwendet, und daher die Umgebung und Terminologie von R verstehen müssen.

Es gibt einige integrierte Entwicklungstools für die Sprache R, die Sie verwenden können. Ein bekanntes Projekt trägt den Namen RStudio. Ich ziehe im Allgemeinen die systemeigene R-Konsole vor. Verschiedene Microsoft-Produkte setzen R ein. Beispielsweise kann der Microsoft Azure Machine Learning-Dienst R-Skripts verwenden, und ich vermute, dass irgendwann einmal Unterstützung für R zu Visual Studio hinzugefügt wird.


Dr. James McCaffrey ist in Redmond (Washington) für Microsoft Research tätig. Er hat an verschiedenen Microsoft-Produkten mitgearbeitet, unter anderem an Internet Explorer und Bing. Dr. McCaffrey erreichen Sie unter jammc@microsoft.com.

Unser Dank gilt den folgenden technischen Experten von Microsoft Research für die Durchsicht dieses Artikels: Dan Liebling und Robin Reynolds-Haertle
Robin Reynolds-Haertle erstellt die Dokumentation für plattformübergreifende Entwicklung mit Visual Studio. Ihre aktuellen Hauptinteressen liegen bei .NET Core, C#, Swift und R