Benötigen Sie einen Algorithmus, um die relative Wahrscheinlichkeit zu berechnen, dass Daten aus der Normal- oder der Lognormalverteilung stammen

13

Angenommen, Sie haben eine Reihe von Werten, und Sie möchten wissen, ob es wahrscheinlicher ist, dass sie aus einer Gaußschen (Normal-) Verteilung oder aus einer logarithmischen Normalverteilung entnommen wurden.

Idealerweise wissen Sie etwas über die Grundgesamtheit oder über die Ursachen von experimentellen Fehlern und hätten daher zusätzliche Informationen, die für die Beantwortung der Frage hilfreich sind. Angenommen, wir haben nur eine Reihe von Zahlen und keine weiteren Informationen. Was ist wahrscheinlicher: Stichproben aus einem Gaußschen oder Stichproben aus einer logarithmischen Normalverteilung? Wie viel wahrscheinlicher? Was ich mir erhoffe, ist ein Algorithmus, mit dem ich zwischen den beiden Modellen wählen und hoffentlich die relative Wahrscheinlichkeit für jedes Modell quantifizieren kann.

Harvey Motulsky
quelle
1
Es könnte eine lustige Übung sein, die Verteilung über Verteilungen in der Natur / veröffentlichten Literatur zu charakterisieren. Andererseits wird es nie mehr als eine lustige Übung sein. Für eine ernsthafte Behandlung können Sie entweder nach einer Theorie suchen, die Ihre Wahl rechtfertigt, oder bei ausreichenden Daten die Anpassungsgüte jeder Kandidatenverteilung visualisieren und testen.
JohnRos
3
Wenn es darum geht, aus Erfahrung zu verallgemeinern, würde ich sagen, dass positiv verzerrte Verteilungen der häufigste Typ sind, insbesondere für Antwortvariablen, die von zentralem Interesse sind, und dass Lognormale häufiger sind als Normalen. Ein Band von 1962, den der Wissenschaftler spekuliert, der vom berühmten Statistiker IJ Good herausgegeben wurde, enthielt ein anonymes Stück "Bloggins 'Arbeitsregeln" mit der Behauptung "Die logarithmische Normalverteilung ist normaler als normal". (Einige der anderen Regeln sind stark statistisch.)
Nick Cox
Ich scheine deine Frage anders zu interpretieren als JohnRos und anxoestevez. Für mich klingt Ihre Frage nach einer einfachen Modellauswahl , dh nach der Berechnung von , wobei entweder die Normalverteilung oder die logarithmische Normalverteilung ist und Ihre Daten sind. Wenn die Modellauswahl nicht das ist, wonach Sie suchen, können Sie dies klären? P(MD)MD
Lucas
@ Lucas Ich denke, deine Interpretation unterscheidet sich nicht so sehr von meiner. In beiden Fällen müssen Sie Apriori- Annahmen treffen .
anxoestevez
2
Warum nicht einfach das verallgemeinerte Wahrscheinlichkeitsverhältnis berechnen und den Benutzer warnen, wenn es das Protokoll-Normal bevorzugt?
Scortchi

Antworten:

7

Sie können den Verteilungstyp am besten schätzen, indem Sie jede (normale oder logarithmische) Verteilung nach maximaler Wahrscheinlichkeit an die Daten anpassen und dann die logarithmische Wahrscheinlichkeit unter jedem Modell vergleichen - das Modell mit der höchsten logarithmischen Wahrscheinlichkeit ist die beste Anpassung. Zum Beispiel in R:

# log likelihood of the data given the parameters (par) for 
# a normal or lognormal distribution
logl <- function(par, x, lognorm=F) {
    if(par[2]<0) { return(-Inf) }
    ifelse(lognorm,
    sum(dlnorm(x,par[1],par[2],log=T)),
    sum(dnorm(x,par[1],par[2],log=T))
    )
}

# estimate parameters of distribution of x by ML 
ml <- function(par, x, ...) {
    optim(par, logl, control=list(fnscale=-1), x=x, ...)
}

# best guess for distribution-type
# use mean,sd of x for starting parameters in ML fit of normal
# use mean,sd of log(x) for starting parameters in ML fit of lognormal
# return name of distribution type with highest log ML
best <- function(x) {
    logl_norm <- ml(c(mean(x), sd(x)), x)$value
        logl_lognorm <- ml(c(mean(log(x)), sd(log(x))), x, lognorm=T)$value
    c("Normal","Lognormal")[which.max(c(logl_norm, logl_lognorm))]
}

Generieren Sie nun Zahlen aus einer Normalverteilung und passen Sie eine Normalverteilung nach ML an:

set.seed(1)
x = rnorm(100, 10, 2)
ml(c(10,2), x)

Erzeugt:

$par
[1] 10.218083  1.787379

$value
[1] -199.9697
...

Vergleichen Sie die Log-Wahrscheinlichkeit für ML-Anpassungen von Normal- und Log-Normalverteilungen:

ml(c(10,2), x)$value # -199.9697
    ml(c(2,0.2), x, lognorm=T)$value # -203.1891
best(x) # Normal

Versuchen Sie es mit einer lognormalen Verteilung:

best(rlnorm(100, 2.6, 0.2)) # lognormal

Die Zuordnung ist abhängig von n, mean und sd nicht perfekt:

> table(replicate(1000, best(rnorm(500, 10, 2))))

Lognormal    Normal 
        6       994 
> table(replicate(1000, best(rlnorm(500, 2.6, 0.2))))

Lognormal    Normal 
      999         1 
waferthin
quelle
1
Sie müssen die Maximum-Likelihood-Parameterschätzungen weder für die normale noch für die logarithmische Normalität numerisch finden (obwohl dies zeigt, wie Sie die Idee auf den Vergleich anderer Verteilungen verallgemeinern würden). Ansonsten sehr vernünftiger Ansatz.
Scortchi
Ich habe R oder das Konzept der maximalen Wahrscheinlichkeit kaum verwendet, daher hier eine grundlegende Frage. Ich weiß, dass wir die AIC (oder BIC) nicht vergleichen können, wenn wir eine Normalverteilung auf die Daten und nicht auf die Protokolle der Daten anwenden, da die AIC oder BIC nicht vergleichbar wären. Man muss zwei Modelle an einen Datensatz anpassen (ohne Transformationen, ohne Ausreißerausschlüsse usw.), und die Transformation der Daten ändert AIC oder BIC, unabhängig davon, ob der Vergleich falsch ist. Was ist mit ML? Ist dieser Vergleich legitim?
Harvey Motulsky
Wir finden die am besten passenden Normal- und Lognormalverteilungen zu den Daten und berechnen dann die Wahrscheinlichkeit der Beobachtung der Daten unter der Annahme, dass sie aus diesen Verteilungen stammen (die Wahrscheinlichkeit oder p(X|\theta)). Wir transformieren die Daten nicht. Wir drucken die Verteilung aus, für die die Wahrscheinlichkeit, die Daten zu beobachten, am höchsten ist. Dieser Ansatz ist legitim, hat aber den Nachteil, dass wir die Wahrscheinlichkeit des Modells bei gegebenen Daten nicht ableiten p(M|X), dh die Wahrscheinlichkeit, dass die Daten aus einer Normal-gegen-Lognormal-Verteilung stammen (z. B. p (normal) = 0,1, p (lognormal) = 0.9) im Gegensatz zum Bayes'schen Ansatz.
Wafer
1
@ Harvey Richtig, aber irrelevant - Sie haben gefragt, ob Sie normale oder logarithmische Verteilungen an die gleichen Daten anpassen möchten. Da die Anzahl der freien Parameter für beide Modelle gleich ist, reduziert sich der Vergleich von AICs oder BICs auf den Vergleich von Log-Wahrscheinlichkeiten.
Scortchi
@wannymahoots Jeder vernünftige Vorgänger für einen Bayes'schen Ansatz in diesem Zusammenhang - basierend auf der Schätzung der relativen Wahrscheinlichkeiten, mit denen ein Softwarebenutzer versucht, normale oder logarithmisch normale Daten abzugleichen - wird so wenig aussagekräftig sein, dass ein Ansatz ähnliche Ergebnisse liefert basierend nur auf der Wahrscheinlichkeit.
Scortchi
10

M{Normal,Log-normal}X={x1,...,xN}

P(MX)P(XM)P(M).

Der schwierige Teil ist, die marginale Wahrscheinlichkeit zu bekommen ,

P(XM)=P(Xθ,M)P(θM)dθ.

p(θM)XY={logx1,...,logxNYX,

P(XM=Log-Normal)=P(YM=Normal)i|1xi|.

P(θM)P(σ2,μM=Normal)P(M)

Beispiel:

P(μ,σ2M=Normal)m0=0,v0=20,a0=1,b0=100

enter image description here

Nach Murphy (2007) (Gleichung 203) ist die marginale Wahrscheinlichkeit der Normalverteilung dann gegeben durch

P(XM=Normal)=|vN|12|v0|12b0a0bnaNΓ(aN)Γ(a0)1πN/22N

aN,bN,vNP(μ,σ2X,M=Normal)

vN=1/(v01+N),mN=(v01m0+ixi)/vN,aN=a0+N2,bN=b0+12(v01m02vN1mN2+ixi2).

Ich verwende die gleichen Hyperparameter für die Log-Normalverteilung,

P(XM=Log-normal)=P({logx1,...,logxN}M=Normal)i|1xi|.

0.1P(M=Log-normal)=0.1

enter image description here

der posterior verhält sich so:

enter image description here

N

Bei der Implementierung der Gleichungen wäre es eine gute Idee, mit logarithmischen Dichten anstelle von Dichten zu arbeiten. Aber sonst sollte es ziemlich einfach sein. Hier ist der Code, mit dem ich die Zeichnungen erstellt habe:

https://gist.github.com/lucastheis/6094631

Lucas
quelle
4

Es hört sich so an, als ob Sie nach etwas sehr Pragmatischem suchen, um Analysten zu helfen, die wahrscheinlich keine professionellen Statistiker sind, und etwas benötigen, das sie zu Standarderkundungstechniken wie dem Betrachten von qq-Diagrammen, Dichtediagrammen usw. auffordert.

Führen Sie in diesem Fall einfach einen Normalitätstest (Shapiro-Wilk oder was auch immer) für die Originaldaten und einen für die logarithmisch transformierten Daten durch. Wenn der zweite p-Wert höher ist, aktivieren Sie ein Flag, damit der Analyst die Verwendung einer logarithmischen Transformation in Betracht zieht ? Als Bonus können Sie eine 2 x 2-Grafik des Dichteliniendiagramms und des qqnorm-Diagramms der Rohdaten und der transformierten Daten ausspucken.

Dies wird Ihre Frage nach der relativen Wahrscheinlichkeit technisch nicht beantworten, aber ich frage mich, ob es alles ist, was Sie brauchen.

Peter Ellis
quelle
Klug. Vielleicht ist das genug und vermeidet die Notwendigkeit, Wahrscheinlichkeitsberechnungen zu erklären ... Danke.
Harvey Motulsky