Afsnit 2.5: Poissonmodellen

Jeg slutter dette kapitel af med at sige noget om poissonfordelingen, som I også kender fra sandsynlighedsdelen af calculus. Jeg vil beskrive, hvordan poissonfordelingen kommer frem som en model for tilfældige ankomster. Dernæst vil jeg estimere parameteren i poissonmodellen og lave et konfidensinterval for denne.

Statistisk Model 2.5.1. (Poissonmodellen)

Lad $X$ være et stokastisk antal, som kan antage alle mulige heltallige værdier større end eller lig med nul. Poissonmodellen med parameter $\lambda$ til beskrivelse af $X$ skrives som

$X\sim\text{poisson}(\lambda),\enspace \lambda\geq 0.$

For poissonmodellen gælder der følgende resultater:

$\begin{aligned} & P(X=x)=\frac{\lambda^x}{x!}e^{-\lambda},\enspace x=0,1,2,\ldots, \\ & E(X)=\lambda,\quad \text{Var}(X)=\lambda. \end{aligned}$ Desuden har vi egenskaben, at hvis $X_1,\ldots,X_k$ er uafhængige og poissonfordelte, så er også summen $X_\bullet=X_1+\cdots+X_k$ poissonfordelt. Den præcise formulering er, at hvis $X_i\sim\text{poisson}(\lambda_i),$ $i=1,\ldots,k,$ og de $k$ stokastiske variable er uafhængige, så er $X_\bullet\sim\text{poisson}(\lambda_\bullet)$ , hvor $X_\bullet$ er summen af $X_i$ -erne og $\lambda_\bullet$ er summen af $\lambda_i$ -erne. Poissonfordelingen bruges til at beskrive fordelingen af et antal, der for eksempel kommer tilfældigt i tid eller tilfældigt over område. Berømte eksempler er antal opkald til en telefoncentral inden for et bestemt tidsinterval (blev studeret af den danske matematiker Erlang), dødsfald ved hestespark i preussiske regimenter, hvor data består af dødsfald per år i 10 regimenter der følges i 20 år (1875‐1894) (billedet nedenfor er fra bogen The art of taming and educating the horse), og fordeling af V1- og V2-bomber over London under anden verdenskrig.

Figur med sandsynligheder

Følgende kode producerer en figur med poissonsandsynligheder. Kør koden, og prøv at ændre på værdien af Lambda i koden (advarsel: benyt aldrig "lambda" som variabelnavn i python). I python beregnes poissonsandsynligheder med st.poisson.pmf.

Kan du se, at "kurveformen" af fordelingen ser ud til at nærme sig en bestemt form, hvis du prøver at gøre $\lambda$ større og større ?

Tilfældige ankomster

Jeg beskriver her en model, der kan forklare, hvorfor mange data kan være poissonfordelte. Modellen skal beskrive situationen med ankomster tilfældigt i tid. Betragt tidsintervallet fra 0 til $t.$ Intervallet deles op i $n$ lige store stykker, hvor $n$ er stor. Uafhængigt af hinanden kaster vi i hvert af de $n$ små intervaller en "mønt", hvor sandsynligheden for krone er $\lambda\frac{t}{n}.$ Hvis mønten viser krone, siger vi, at der er en ankomst i det tilsvarende interval. På denne måde er sandsynligheden for en ankomst i et lille interval proportional med længden af intervallet. Proportionalitetskonstanten er $\lambda,$ som får fortolkningen som rate, det vil sige det forventede antal ankomster per tid. Uafhængigheden mellem de små intervaller modellerer, at ankomsterne er tilfældige i tid.

Ideen er nu at lade $n$ blive større og større, svarende til en finere og finere inddeling af intervallet $[0,t]$ i små intervaller. Lad $X$ være den stokastiske variabel, der angiver antal ankomster i intervallet $[0,t].$ Når $n$ er valgt, står vi i en binomialsituation og har $X\sim\text{binom}(n,\lambda t/n).$ Lader vi $n$ blive større og større, kan man matematisk vise, at

$P(X=x)=\binom{n}{x}(\lambda t/n)^x(1-\lambda t/n)^{n-x} \rightarrow \frac{(t\lambda)^x}{x!}e^{-t\lambda}.$ Her står, at i grænsen hvor $n$ bliver meget stor, vil $X$ være poissonfordelt med parameter $t\lambda,$ $X\sim\text{poisson}(t\lambda).$

Man udtrykker ovenstående argument på den måde, at en binomialfordeling $\text{binom}(n,p)$ med $n$ stor og $p$ lille ligner en poissonfordeling med parameter $np.$ Nedenstående kommandovindue laver en figur til illustration af dette. I figuren er poissonsandsynlighederne angivet ved søjler og binomialsandsynlighederne er angivet ved røde plustegn. Kør koden, og prøv dernæst at sætte $n$ til 100. Kan du se hvad der sker ? Prøv også en anden værdi af Lambda.

Test dig selv: valg af statistikmodel

Denne lille testopgave går ud på at træne jer i at kunne vurdere, om data skal beskrives med en binomialmodel eller med en poissonmodel.

Quiz

I løbet af en måned er der indsamlet 20 vandprøver fra forskellige vandboringer, og 5 af disse har pesticidrester over det tilladte niveau. Klik på de rigtige udsagn nedenfor, hvor Pest står for den bagvedliggende stokastiske variabel.

De 5 prøver kommer tilfældigt over tid i løbet af en måned, hvorfor det er en observation fra en poissonfordeling, $\text{Pest}\sim\text{poisson}(\lambda t),$ hvor $t$ er en måned.

Hver af de 20 prøver kan have pesticidrester enten over eller under grænsen, hvorfor 5 er en observation fra en binomialfordeling, $\text{Pest}\sim\text{binom}(20,p).$

Prøv også at forklare, hvorfor det forkerte valg er forkert.

Svar: Forkerte statistikmodel

Der er ikke tale om, at der kommer prøver med højt indhold af pesticidrester tilfældigt over tid. Tidsaspektet i eksperimentet har ikke noget at gøre med indeholdet af pesticidrester, men siger blot noget om, hvornår man vælger at tage ud og indsamle vandprøven.

2.5.1 Estimation i poissonmodellen

Betragt uafhængige poissonfordelte variable $X_1,\ldots,X_n,$ $X_i\sim\text{poisson}(t\lambda),$ hvor $t$ er kendt og $\lambda>0$ er en ukendt parameter, som vi ønsker at estimere. Bemærk her, at parameteren i poissonfordelingen er $t\lambda$ , hvor vi tænker på $t$ som et tidsrum, og $\lambda$ fortolkes som den forventede rate per tid. På grund af uafhængigheden bliver likelihoodfunktionen et produkt af poissonsandsynligheder

$L(\lambda)=\prod_{i=1}^n \frac{(t\lambda)^{x_i}}{x_i!}e^{-t\lambda}= \Big(\frac{1}{\prod_ix_i!} \Big) (t\lambda)^{x_1+\cdots+x_n}e^{-nt\lambda}.$ For at finde den værdi af $\lambda$ som giver maksimum af likelihoodfunktionen, tager vi logaritmen, differentierer med hensyn til $\lambda$ og sætter den afledede lig med nul. Dette giver

$\frac{x_1+\cdots+x_n}{\lambda}-nt=0\quad\text{eller}\quad \hat\lambda=\frac{\bar x}{t},\quad \bar x=\frac{x_1+\cdots+x_n}{n}.$ Ligesom i binomialmodellen giver dette skøn god mening: $\lambda$ er den forventede rate per tid, og $\hat\lambda$ er den observerede rate per tid.

Eksempel 2.5.2. (Rutherford og Geiger)

I en berømt artikel fra 1910 giver Rutherford og Geiger resultaterne fra 2608 målinger af antal henfald fra en poloniumsmasse i 72 sekunder. Vi kan beskrive situationen med $X_i\sim\text{poisson}(72\cdot\lambda),$ $i=1,\ldots,2608.$ I denne situation bliver skønnet over raten per sekund $\hat\lambda=10097/(72\cdot 2608)=0.05377,$ idet summen af de 2608 målinger er 10097. Et histogram med de 2608 målinger er vist i den følgende figur.

2.5.2 Konfidensinterval i poissonmodellen

I poissonmodellen $X\sim\text{poisson}(t\lambda)$ er det ikke muligt at lave et konfidensinterval for $\lambda,$ som opfylder betingelsen i Definition 2.2.1 eksakt. I stedet benytter man følgende approksimative konfidensinterval.

Resultat 2.5.3. (Konfidensinterval i poissonmodellen)

For modellen $X\sim\text{poisson}(t\lambda)$ er et approksimativt 95%-konfidensinterval for $\lambda$ baseret på observationen $x$ givet ved

$\Big[\frac{1}{t}\Big(x+\frac{u^2}{2}-u\sqrt{x+\frac{u^2}{4}}\Big),\, \frac{1}{t}\Big(x+u^2/2+u\sqrt{x+u^2/4}\Big)\Big],\quad u=1.96.$

Hvis værdien af $u$ ændres til 1.645, får man i stedet et approksimativt 90%-konfidensinterval. Ligesom i binomialmodellen bygger konfidensintervallet her på den centrale grænseværdisætning.

Eksempel 2.5.4. (Rutherford og Geiger)

I Eksempel 2.5.2 omkring henfald af polonium er der 2608 uafhængige målinger: $X_i\sim\text{poisson}(72\cdot\lambda),$ $i=1,\ldots,2608.$ Som nævnt ovenfor er summen af $X_i$ -erne også poissonfordelt. Vi har derfor, at $X=\sum_i X_i\sim\text{poisson}(72\cdot 2608\cdot\lambda).$ Den målte værdi af $X$ er $x=10097.$ I kodevinduet nedenfor udregnes det approksimative 95%-konfidensinterval.

2.5.5 Beregning i python af konfidensinterval for rate

Den følgende kode kan bruges generelt til beregning af konfidensinterval for en rate, hvis man selv indskriver det observerede antal ( $x$ ) og tidsparameteren ( $t$ ).

Foregående Næste