Verwendung des Chi-Quadrat-Tests, um festzustellen, ob die Daten der Poisson-Verteilung folgen

Die folgende Abbildung (Abbildung 1 aus S. 646 dieses Papiers ) vergleicht die beobachteten Werte mit den erwarteten Werten unter der Poisson-Verteilung. Anschließend wird ein Chi-Quadrat-Test durchgeführt, um festzustellen, ob die beobachteten Werte von den erwarteten Werten unter der Poisson-Verteilung abweichen.

Wie ist es mit R möglich, erwartete Werte unter Poisson-Verteilung zu generieren und beobachtete Werte mit einem Chi-Quadrat-Test zu vergleichen?

BEARBEITEN:

Hier ist mein Versuch, das zu tun, was sie in Papierform getan haben. Ich möchte wissen, ob sich die beobachtete Verteilung von variablevon einer Poisson-Verteilung unterscheidet. Ich möchte auch wissen, ob das, was ich unten getan habe, das gleiche Verfahren ist wie das, was sie in Papierform getan haben. Da der P-Wert> 0,05 ist, habe ich unten festgestellt, dass die Verteilung von variableeiner Poisson-Verteilung folgt - könnte jemand dies bestätigen?

df <- data.frame(variable = 0:5, frequency = c(20, 10, 5, 3, 2, 1))

# estimate lambda
mean_df_variable <- mean(df$variable)

# calculate expected values if df$frequency follows a poisson distribution
library(plyr)
expected <- laply(0:5, function(x) dpois(x=x, lambda=mean_df_variable, log = FALSE))

# calculate actual distribution of df$frequency
observed <- df$frequency/sum(df$frequency)

# does distribution of df$frequency differ from a poisson distribution? Apparently 
#   not because P-value is > 0.05
chisq.test(expected, observed)

Die Art und Weise, wie Sie den Chi-Quadrat-Test durchgeführt haben, ist nicht korrekt. Es gibt mehrere Probleme. Zunächst sieht Ihr Datenrahmen folgendermaßen aus:

  variable frequency
1        0        20
2        1        10
3        2         5
4        3         3
5        4         2
6        5         1

Wenn du also rennst mean(df$variable), bekommst du 2.5, was nur der Mittelwert von ist 0:5. Das heißt, es ist ungewichtet. Erstellen Sie stattdessen Ihre Variable wie folgt:

x = rep(0:5, times=c(20, 10, 5, 3, 2, 1))
table(x)
# x
#  0  1  2  3  4  5 
# 20 10  5  3  2  1
mean(x)
# [1] 1.02439

Der table()Aufruf zeigt, dass der Code uns das gibt, was wir wollten, und mean()schätzt Lambda daher korrekt.

Als nächstes gehen Ihre geschätzten Wahrscheinlichkeiten nur auf 5, aber die Poisson-Verteilung geht auf unendlich. Sie müssen also die Wahrscheinlichkeiten der Werte berücksichtigen, die Sie nicht in Ihrem Datensatz haben. Dies ist nicht schwer zu tun, Sie berechnen nur die Ergänzung:

probs = dpois(0:5, lambda=mean(x))
probs
# [1] 0.359015310 0.367771781 0.188370912 0.064321775 0.016472650 0.003374884
comp = 1-sum(probs)
# [1] 0.0006726867

Schließlich sind in und R's chisq.test()die Argumente x=und y=nicht genau für die erwarteten und beobachteten Werte in der Art und Weise, wie Sie dies einrichten. Zum einen nennen Sie "erwartet" tatsächlich Wahrscheinlichkeiten (dh die Ausgabe von dpois()). Um diese erwarteten Werte zu erhalten, müssten Sie diese Wahrscheinlichkeiten (und das Kompliment) mit der Gesamtzahl multiplizieren. Aber selbst dann würden Sie diese nicht verwenden y=. Auf jeden Fall müssen Sie das eigentlich nicht tun, Sie können dem p=Argument einfach die Wahrscheinlichkeiten zuweisen . Darüber hinaus müssen Sie 0Ihrem beobachteten Wertevektor einen hinzufügen , um alle möglichen Werte darzustellen, die in Ihrem Datensatz nicht angezeigt werden:

chisq.test(x=c(20, 10, 5, 3, 2, 1, 0), p=c(probs, comp))

#  Chi-squared test for given probabilities
# 
# data:  c(20, 10, 5, 3, 2, 1, 0)
# X-squared = 12.6058, df = 6, p-value = 0.04974
# 
# Warning message:
#   In chisq.test(x = c(20, 10, 5, 3, 2, 1, 0), p = c(probs, comp)) :
#   Chi-squared approximation may be incorrect

Die Warnmeldung deutet darauf hin, dass wir möglicherweise lieber simulieren möchten, also versuchen wir es erneut:

chisq.test(x=c(20, 10, 5, 3, 2, 1, 0), p=c(probs, comp), simulate.p.value=TRUE)

# Chi-squared test for given probabilities with simulated p-value 
#   (based on 2000 replicates)
# 
# data:  c(20, 10, 5, 3, 2, 1, 0)
# X-squared = 12.6058, df = NA, p-value = 0.07046

Dies ist vermutlich ein genauerer p-Wert, wirft jedoch die Frage auf, wie er interpretiert werden sollte. Sie fragen: "Da der P-Wert> 0,05 ist, habe ich unten festgestellt, dass die Verteilung der Variablen einer Poisson-Verteilung folgt - könnte jemand dies bestätigen?" Bei Verwendung des richtigen Ansatzes stellen wir fest, dass der erste p-Wert nur <0,05 war, der zweite (simulierte) p-Wert jedoch nur> 0,05. Obwohl der letztere p-Wert genauer ist, würde ich nicht schnell zu dem Schluss kommen, dass die Daten aus einer Poisson-Verteilung stammen. Hier sind einige Fakten zu beachten:

• Wie im Titel einer Arbeit von Gelman und Stern angegeben, ist der Unterschied zwischen "signifikant" und "nicht signifikant" selbst statistisch nicht signifikant .
• Reale Daten stammen nie aus einer der idealisierten Verteilungen, eine Tatsache, auf die @Glen_b kürzlich hier gut hingewiesen hat: Wie sieht dieses Histogramm aus?

  Reale Daten tendieren dazu, nicht den einfachen Verteilungsformen der üblichen Ein-, Zwei- oder Drei-Parameter-Verteilungen zu folgen. Reale Verteilungen ähneln eher heterogenen Gemischen. Einfache Verteilungsformen sind bequeme Fiktionen (genauer gesagt Modelle) - sie nähern sich der Realität auf eine Weise an, die das Arbeiten erleichtert.

• Sie können die Tatsache eines nicht signifikanten Ergebnisses nicht verwenden, um die Nullhypothese zu bestätigen, wie ich hier erkläre: Warum sagen Statistiker, dass ein nicht signifikantes Ergebnis bedeutet, dass Sie die Nullhypothese nicht ablehnen können, anstatt die Nullhypothese zu akzeptieren?

Wenn ich verstanden habe, was Sie meinten, sollten Sie:

1. Schätzen Sie den Parameter der Poisson-Verteilung für Ihre Daten, vorausgesetzt, es handelt sich beispielsweise um eine Poisson-Verteilung

lambdaEst = mean(x)

2. $0,1, 2, ...$

probTheo0 = dpois(x = 0, lambda = lambdaEst, log = FALSE)

3. Vergleichen Sie dann die tatsächlichen mit den theoretischen Wahrscheinlichkeiten über einen Chi-Quadrat-Test nach diesem Ansatz. ChiSquare Test CV-Lösung