Mennyit teszteljünk? 1. rész

Napjainkban sok szó esik arról, hogy hány tesztet kellene elvégezni ahhoz, hogy kellő pontossággal ismerni lehessen a fertőzöttek vagy fertőzésen átesettek számát egy bizonyos betegség esetén. Ez számos matematikai kérdést felvet. Ebben az írásban arról szólunk, hogy ideális feltételek mellett a fertőzésen átesettek országos arányát hány teszt elvégzésével lehet biztonságosan megállapítani.

Egyszerű számpéldával a következőképpen lehet a kérdést megvilágítani. Tegyük fel, hogy a társadalom $10\%$ -a esett át a fertőzésen, de ennek kiderítéséhez természetesen nem szeretnénk 10 millió emberen a tesztet elvégezni, csak kisebb mintán. Például, ha 100 embertől veszünk mintát, akkor, ha szerencsénk van, éppen 10 fertőzésen átesett egyént találunk, ami éppen az országos átlagot adná, de könnyen lehet, hogy 100-ból 8-at vagy 11-et kapunk, és tévesen becsüljük meg a valódi arányt. Ekkor tesztelhetünk újabb száz embert, majd ismét újabbat. A kérdés az, hogy hányszor kell ezt elvégezni, és milyen biztonsággal tudjuk a példa elején említett 10 százalékot megkapni. Ezt a bevezető valószínűségszámítás kurzusokon (közgazdasági, mérnöki szakokon is) elsajátított matematikai eszköztár segítségével magunk is kiszámolhatjuk, ezért örömmel osztjuk meg az Olvasóval ennek egy lehetséges módszerét. A végeredményt még nem áruljuk el, a kíváncsiak előrelapozhatnak és megleshetik. Aki pedig sorban haladva olvassa a cikket, azt elkalauzoljuk olyan témákhoz, mint a Stirling-formula, Gauss-görbe és centrális határeloszlás-tétel, sőt még az is kiderül, hogy nem volt felesleges megismerni analízisből az integrálközelítő összeg fogalmát és megtanulni a helyettesítéses integrálást.

1. A mintavételezés modellje: a binomális eloszlás

Induljon hát a matematikaóra. Ahhoz, hogy egy adott pontossághoz és bizonyossághoz tartozó szükséges tesztek számát meghatározzuk, először azt a kérdést értjük meg alaposabban, hogy adott számú teszt és ismert védettségi arány esetén mennyi a valószínűsége, hogy például legfeljebb $r=5\%$ -ot téved a becslésünk az igazi érték arányában. A számolások során azt feltételezzük, hogy a tesztek egymástól függetlenek, és tökéletes megbízhatósággal mutatják ki a fertőzést, nincs tévedés egyik irányban sem. A kérdés az egyszerű pénzfeldobáshoz (Bernoulli trials) hasonló, az egyének $p$ valószínűséggel átestek a fertőzésen, őket védetteknek fogjuk nevezni, $1-p$ valószínűséggel pedig nem. Ennek a $p$ -nek az értékre leszünk kíváncsiak adott pontossággal. Legyen $X$ azon emberek száma $n$ vizsgált személy közül, akiknél kimutatta a teszt a védettséget. Ennek értéke véletlen, azaz $X$ egy valószínűségi változó. Ekkor a $p$ -t becsülni az $X/n$ -nel, vagyis a védettek relatív gyakoriságával tudjuk a legegyszerűbben és sok szempontból a legjobban (egy ezzel kapcsolatos kérdésre a cikk második részében fogunk kitérni). A becslésre azt a feltételt fogalmaztuk meg, hogy a relatív hibája legfeljebb $r=5\%$ legyen:

$\displaystyle \bigg\vert\frac Xn-p\bigg\vert\leq 0{,}05p.$

(1)

Mivel a védettek száma, azaz $X$ értéke 0 és $n$ között tetszőleges lehet, olyan feltételt nem fogunk tudni mondani, amire ez biztosan teljesül, bárhogyan is alakul a mintavételezés. Azt azonban már kitűzhetjük célnak, hogy ez minél nagyobb, például legalább $95\%$ valószínűséggel teljesüljön.

Ehhez pedig először a fenti esemény valószínűségét szeretnénk tehát kiszámítani. Ha $n$ egyént tesztelünk, akkor a középiskolában tanult képlet alapján annak valószínűsége, hogy $k$ védettet találunk:

$\displaystyle p_k= \mathbb{P}(X=k)=\binom nk p^k (1-p)^{n-k} \qquad (k=0, 1, 2, \ldots, n).$

(2)

Ez a jól ismert binomiális eloszlás, ami azt is megadja, hogy ha egy urnában a golyók $p$ -ed része piros, akkor mennyi annak valószínűsége, hogy $n$ -szer visszatevéssel húzva pontosan $k$ pirosat kapunk.

1. ábra. A binomiális eloszlás $n=1000$ renddel és $p=0{,}1$ paraméterrel, valamint a közelítő függvény (4) alapján

Már az egyetemi bevezető matematika órára átlépve, tanultuk, hogy ennek várható értéke $np$ , azaz átlagosan ennyi védettet fogunk találni az $n$ tesztelt egyén között. Ha például $n=1000$ tesztelést végzünk, és az emberek 10 százaléka védett, azaz $p=0{,}1$ , akkor várhatóan 100 védettet fogunk találni. Valójában persze lehet, hogy 95-öt vagy 110-et. Ha 95 védettet találunk, akkor a becslésünk $X/n=95/1000=0{,}095$ lesz, ezzel az igazi $p$ -hez viszonyítva legfeljebb $5\%$ -t tévedtünk, azaz az $r=5\%$ relatív hibahatáron belül maradunk. Ha viszont $X=110$ , akkor a becslésünk 0,11, ez már $10\%$ -os hibát jelentene. Vagyis az (1) egyenlet alapján annak a valószínűségét szeretnénk kiszámítani, hogy a $95$ -től $105$ -ig terjedő intervallumba esik a védettek száma. A (2) képletből ez a valószínűség egyszerűen a

$\displaystyle p_{95} + p_{96}+ \ldots + p_{105}$

összeg. Az 1. ábrán láthatjuk a binomiális eloszlást $k=80$ -tól $k=120$ -ig, $n=1000$ és $p=0{,}1$ esetén. Az oszlopok szélessége $1$ , a magasságokat a fenti $p_k$ -k adják meg, így a sötétebb rész területét kell meghatároznunk.

Az alábbi octave kód abban lesz segítségünkre, hogy ezt numerikusan kiszámítsuk. Más programokban, például Matlabban vagy R-ben is egy hasonlóképpen megírt algoritmust használhatnánk, illetve ezen programoknak a binomiális eloszlásra vonatkozó beépített függvényeit is (például binocdf(105, 1000, 0.1)-binocdf(94, 1000, 0.1) az octave-ban). Ez a hosszabb kód viszont a számítás menetét is illusztrálja. Az mv vektor azokat a lehetséges értékeket tartalmazza, melyek a megengedett tartományba esnek. A pk függvény adott $n, p$ számok és $k$ , pozitív egészekből álló vektor esetén egy olyan vektort ad vissza, aminek a $k.$ eleme éppen a fent definiált $p_k$ . A logaritmusra való áttérés oka, hogy nagy $n$ esetén a binomiális együttható értéke is nagyon nagy, $p^k$ pedig nagyon kicsi lehet, és a szorzásnál megbízhatóbb a logaritmusok kiszámítása és összeadása, így kisebb a számítógép kerekítéseiből adódó hiba. Végül a keresett tartományba eső valószínűségeknek, azaz binmv elemeinek az összege adja meg a kérdéses valószínűséget.

n=1000; p=0.1; hiba=0.05;

mv=[round(n*p*(1-hiba)): round(n*p*(1+hiba))];

binmv=bineloszl(n,p,mv);

function pk=bineloszl(n,p,k)

bnc=[];

for j=1:length(k)

mj=m(j);

bnc=[bnc bincoeff(n,mk)];

end

lpk=log(bnc)+k*log(p)+(n-k)*log(1-p);

pk=exp(lpk);

endfunction

sum(binmv)

Ebből azt kapjuk, hogy $n=1000$ és $p=0{,}1$ esetén annak valószínűsége, hogy a becslésünk relatív hibája legfeljebb $r=5\%$ , vagyis az (1) képletben felírt esemény valószínűsége:

$\displaystyle \mathbb{P}\bigg(\bigg\vert\frac Xn-p\bigg\vert\leq 0{,}05p\bigg)=... ...\vert X-np\vert\leq 0{,}05\cdot np\big)=\mathbb{P}(95\leq X\leq 105)=43{,}79\%.$

(3)

2. A binomiális eloszlás tagjainak közelítése a haranggörbével

Ezzel tehát adott mintaelemszám esetén már tudunk számolni. Ugyanakkor az imént kapott érték azt jelenti, hogy több mint $50\%$ valószínűséggel tévedünk a megengedett $5\%$ relatív hibánál többet. Ennél jobb eredményt szeretnénk elérni, vagyis több mint ezer tesztre van szükség. A fenti közvetlen számolásból azonban nem látszik, hogy pontosan hány darabra, ráadásul nagyobb $n$ esetén a közelítésekből túlságosan sok hiba adódhat (például ha $n=3000$ és $k=40$ , akkor a binomiális együttható már $91$ -jegyű, ha $k=400$ , akkor több mint $1200$ számjegyből áll). Ennél egyszerűbben és megbízhatóbban számolhatunk, ha a fenti módszer helyett, ahogy az 1. ábra is sugallja, a kapott valószínűséget egy integrálközelítő összegnek tekintjük.

Így ahhoz a kérdéshez jutunk, hogy mi is az a függvény, amit ezek a $p_k$ valószínűségek közelítenek. Mivel a binomiális együtthatókban faktoriálisok szerepelnek, a közelítésre az ezek aszimptotikus viselkedését leíró Stirling-formula ad lehetőséget:

$\displaystyle \lim_{n\rightarrow\infty} \frac{n!}{\big(\frac ne\big)^n \sqrt{2\pi n}}=1.$

Innen elindulva számos további lépésen keresztül juthatunk el a de Moivre–Laplace-tétel lokális alakjáig, mely szerint

$\displaystyle p_k\approx \frac{1}{\sqrt{np(1-p)}}\cdot\frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(k-np)^2}{2np(1-p)}\bigg),$

(4)

ha a $k-n$ távolság nem több $n^{2/3}$ -nál, $n$ megfelelően nagy, illetve $p$ és $1-p$ közül egyik sem túlságosan kicsi. A tétel ennél valójában sokkal erősebbet is állít, és pontosabban is megfogalmazható [4], [8]. Hasonló állítást de Moivre már 1730-ban publikált a $p=1/2$ esetben [7], Laplace pedig ezt általánosította a XIX. század első felében más $p$ értékekre [6]. Ezzel tehát kapunk egy képletet az 1. ábrán látható tengelyesen szimmetrikus közelítő függvényre. Az ehhez hasonló (ebből eltolással vagy nyújtással megkapható) függvényeket szokták haranggörbének vagy Gauss-görbének nevezni.

Az alábbi octave kód azt számolja, hogy mit kapunk a $\mathbb{P}(95\leq X\leq 105)$ valószínűségre, ha az előző esethez hasonlóan először meghatározzuk a megfelelő értékek halmazát, mindegyikre kiszámítjuk a $p_k$ -nak a most megadott közelítését, majd ezeket összeadjuk.

n=10000; p=0.1; hiba=0.05;

mv=[round(n*p*(1-hiba)): round(n*p*(1+hiba))];

function pkkoz=bineloszlnormkoz(n,p,kv)

kvu=-1*(kv-n*p).2/(2*n*p*(1-p));

pkkoz=exp(kvu)/(sqrt(pi*2*n*p*(1-p)));

endfunction

binkoz=bineloszlnormkoz(n,p,mv);

sum(binkoz)

A számolást elvégezve az eredmény $43{,}81\%$ lesz. A (3) egyenlettel összehasonlítva láthatjuk, hogy ez valóban jó közelítést ad a kérdéses valószínűségre.

3. Összegzés helyett integrálás

Az 1. ábra alapján az volt a célunk, hogy a sötétebb rész területét, azaz a $\mathbb{P}(95\leq X\leq 105)$ valószínűséget a haranggörbe alatti területtel közelítsük. Pontosabban, a (4) egyenlet alapján (figyelve arra, hogy az integrálás is egy $11$ hosszú intervallumon történjen):

$\displaystyle \mathbb{P}(95\leq X\leq 105)$	$\displaystyle =\sum_{k=95}^{105} p_k\approx \sum_{k=95}^{105} \frac{1}{\sqrt{np... ...)}}\cdot \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(k-np)^2}{2np(1-p)}\bigg)\approx$
	$\displaystyle \approx\int_{94{,}5}^{105{,}5} \frac{1}{\sqrt{np(1-p)}}\cdot \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{(y-np)^2}{2np(1-p)}\bigg)dy =$
	$\displaystyle =\int_{94{,}5-np}^{105{,}5-np} \frac{1}{\sqrt{np(1-p)}}\cdot \frac{1}{\sqrt{2\pi}}\exp\bigg(-\frac{z^2}{2np(1-p)}\bigg)dz=$
	$\displaystyle =\int_{\frac{94{,}5-np}{\sqrt{np(1-p)}}}^{\frac{105{,}5-np}{\sqrt{np(1-p)}}} \frac{1}{\sqrt{2\pi}} \exp\bigg(-\frac{x^2}{2}\bigg)dx,$

ahol a $z=y-np$ , majd az $x=\frac{z}{\sqrt{np(1-p)}}$ változócseréket hajtottuk végre, hogy a függvény mindig ugyanaz legyen, csak a határok függjenek a változó paraméterektől.

A következő lépés az lehetne, hogy az integrálra explicit képletet adjunk az $\exp(-x^2/2)$ primitív függvényének segítségével. Ilyen azonban sajnos nincs, be is lehet bizonyítani, hogy az alapvető műveletek segítségével felírható függvények között nem találhatunk megfelelőt. Ezért az alábbi szokásos jelölést vezetjük be arra, hogy mennyi az integrál $t$ -ig:

$\displaystyle \Phi(t)=\int_{-\infty}^t \frac{1}{\sqrt{2\pi}} \exp\bigg(-\frac{x^2}{2}\bigg)dx.$

Ennek a függvénynek, amit a standard normális eloszlás eloszlásfüggvényének neveznek (erre még később visszatérünk), a $t$ helyen felvett értéke a 2. ábrán a $t$ -től balra eső terület nagysága. Így például az is leolvasható, hogy $\Phi(1)=0{,}001+0{,}021+0{,}138+0{,}341\cdot 2=0{,}842$ . A teljes görbe alatti terület pedig $1$ (azt, hogy $t\rightarrow\infty$ esetén $1$ az integrál, polárkoordinátákra való áttéréssel lehet bizonyítani például, de ezt most nem tesszük meg, erről bővebben például itt lehet olvasni: [7.3. szakasz, 3], a normális eloszlásról pedig a [2] szimulációkat tartalmazó jegyzetben is).

2. ábra. Az $\frac{1}{\sqrt{2\pi}}e^{-x^2/2}$ függvény és az alatta lévő területek

A fenti számolás általánosításával oda jutunk, hogy megfelelő feltételeket teljesítő $a_n, b_n$ számsorozatokra:

$\displaystyle \mathbb{P}(a_n \leq X\leq b_n)\approx\Phi\bigg(\frac{b_n-np}{\sqrt{np(1-p)}}\bigg)-\Phi\bigg(\frac{a_n-np}{\sqrt{np(1-p)}}\bigg),$

(5)

hiszen a különbség nem más, mint a két érték közötti integrál. A példában pedig (az octave-ban $\Phi(x)$ a normcdf(x) paranccsal kapható meg):

$\displaystyle \mathbb{P}(95 \leq X\leq 105)\approx\Phi\bigg(\frac{5{,}5}{\sqrt{... ...\Phi\bigg(-\frac{5{,}5}{\sqrt{90}}\bigg)=\Phi(0{,}58)-\Phi(-0{,}58)=43{,}74\%.$

A (3) egyenlettel összehasonlítva látható, hogy továbbra is jó közelítést kaptunk.

4. Kapcsolat a centrális határeloszlás-tétellel és a normális eloszlással

Ezt a közelítő számolást precízzé téve juthatnánk el a de Moivre–Laplace-tétel globális alakjához (vagy annak általánosabb formájához, a centrális határeloszlás-tételhez [2], [3], [4], [8], [9]), és ebből a „megfelelő feltételek” is kiderülnek. A tételből ugyanis az következik, hogy rögzített $u, v$ valós számok és rögzített $0<p<1$ esetén, ha $X_n$ binomiális eloszlású $n$ renddel és $p$ paraméterrel, akkor

$\displaystyle \lim_{n\rightarrow\infty} \mathbb{P}(np+u\sqrt{np(1-p)}\leq X_n\leq np+v\sqrt{np(1-p)})=\Phi(v)-\Phi(u).$

Ha tehát rögzített $u$ és $v$ mellett

$\displaystyle a_n=np+u\sqrt{np(1-p)}$ és $\displaystyle \qquad b_n=np+v\sqrt{np(1-p)},$

akkor a fenti közelítésben $n\rightarrow\infty$ határértéket véve pontos egyenlőséget kapunk. Ebből az is látható, hogy ezzel a módszerrel a várható értéktől való $\sqrt n$ nagyságrendű eltérések valószínűségét tudjuk jól meghatározni (ami éppen a binomiális eloszlás szórásának a nagyságrendje). A centrális határeloszlás-tétel abban általánosabb, hogy nem csak a binomiális eloszlásról, hanem tetszőleges független azonos eloszlású, véges szórású valószínűségi változók összegéről szól, és arról látja be, hogy a várható értéktől való $\sqrt n$ nagyságrendű eltérések limesze a fentihez hasonlóan a standard normális eloszlásfüggvénnyel írhatók le.

A közelítés hibája is felülről korlátozható, erre többféle lehetőség is van. Ennek egy egyszerű, bár nem a legjobb becslést adó alakja, ami a centrális határeloszlás egy még erősebb változatából, a Berry–Esséen-tételből [4] kapható, azt mondja, hogy a fenti közelítésénél (csak az egyik tagot tekintve) legfeljebb $\frac{0{,}48(p^2+(1-p)^2)}{\sqrt {np(1-p)}}$ hibát követünk el. Ebből az mindenképpen fontos, hogy ha $np(1-p)$ közel van nullához, azaz $p$ vagy $1-p$ túl kicsi az $n$ -hez képest, akkor a közelítés elromolhat. Ezért ennek a módszernek az alkalmazásánál például azt a feltételt szokták előírni, hogy $np(1-p)\geq 10$ teljesüljön [9]. Mivel $X$ várható értéke $np$ , ez nem jelent sokkal többet, mint hogy legyen akkora az $n$ , hogy legalább $10-20$ védett embert találjunk.

Végül röviden arról, hogy mi is a standard normális eloszlás, aminek eloszlásfüggvénye $\Phi$ . Ez egy olyan valószínűségi változó, amire az igaz tetszőleges $u<v$ esetén, hogy $\Phi(v)-\Phi(u)$ valószínűséggel esik az $u$ és $v$ számok közé, és ilyen módon $\Phi(v)$ valószínűséggel lesz kisebb $v$ -nél. Ahhoz, hogy egy ilyen $Z$ valószínűségi változót kapjunk, nem kell mást tenni, mint a 2. ábrán látható görbe alatti, $1$ nagyságú területről kiválasztani egy pontot „egyenletesen”, véletlenszerűen, majd $Z$ -nek ennek a pontnak a vízszintes tengelyre eső koordinátáját választani. A centrális határeloszlás-tétel pedig úgy is fogalmazható, hogy az $(X_n-np)/\sqrt{np(1-p)}$ valószínűségi változó „hasonlóan viselkedik” a standard normális eloszláshoz.

5. A szükséges tesztek száma

Az előzőek alapján már közelítő választ tudunk adni az alábbi kérdésre.

Kérdés. Tegyük fel, hogy $n$ embert megvizsgálva mindenki egymástól függetlenül $p$ valószínűséggel védett. Mekkorára válasszuk $n$ -t, hogy annak valószínűsége, hogy az elvégzett tesztek eredménye alapján adott becslésünknek és $p$ -nek a különbsége legfeljebb $h=1\%$ , legalább $q=95\%$ legyen?

Legyen továbbra is $X$ az, hogy az $n$ elvégzett teszt közül hány mutat ki védettséget. A becslésünk továbbra is az $X/n$ relatív gyakoriság lesz. A kérdés feltételét a következőképpen fogalmazhatjuk meg, az (1) egyenletből kiindulva:

$\displaystyle \mathbb{P}\bigg(\bigg\vert\frac Xn-p\bigg\vert\leq 0{,}01\bigg)\geq 0{,}95$ minden $\displaystyle 0<p<1$ -re $\displaystyle .$

A feltétel kicsit más, mint az előző részben, ott ugyanis a relatív hibát adtuk meg ( $p$ -nek hányadrésze lehet az eltérés), most viszont az abszolút hibát (legfeljebb milyen messze lehet a becslés $p$ -től). A relatív hiba esetére később még visszatérünk, ahogy érintőlegesen arra az esetre is, ha esetleg $p$ -ről van valamilyen előzetes információnk. A kérdés egyelőre abból indul ki, hogy olyan módszert keresünk, ami az ismeretlen $p$ valószínűség minden értékére jó, holott természetesen nem biztos, hogy ugyanaz a módszer lenne egyformán hatékony egy ezrelékes vagy tízszázalékos védettségi arány esetén.

A fenti valószínűséget alakítsuk át a korábban látott (5) közelítésnek megfelelően (az $a_n, b_n$ -re vonatkozó „megfelelő feltételeket” egyelőre figyelmen kívül hagyva):

$\begin{equation*}\begin{aligned}\mathbb{P}\bigg(\bigg\vert\frac Xn-p\bigg\vert\l... ...i\bigg(\frac{0{,}01\sqrt n}{\sqrt{p(1-p)}}\bigg)-1. \end{aligned}\end{equation*}$

(6)

Az utolsó lépés abból adódik, hogy a 2. ábrán látható haranggörbe 0-ra való szimmetriája miatt a $-t$ -től balra lévő terület, vagyis $\Phi(-t)$ ugyanaz, mint a $t$ -től jobbra lévő, mivel viszont a teljes görbe alatti terület $1$ , ez utóbbi nem más, mint $1-\Phi(t)$ .

Tehát azt szeretnénk elérni, hogy minden $0<p<1$ -re

$\displaystyle 2\Phi\bigg(\frac{0{,}01\sqrt n}{\sqrt{p(1-p)}}\bigg)-1\geq 0{,}95... ...trightarrow\quad\Phi\bigg(\frac{0{,}01\sqrt n}{\sqrt{p(1-p)}}\bigg)\geq 0{,}975$

(7)

teljesüljön. Vagyis a törtnek olyan értéket kell felvennie, amire igaz, hogy a 2. ábrán a tőle balra eső terület legalább $0{,}975$ . Az ábrán azt látjuk, hogy a $2$ -től balra eső terület éppen 1 – 0,022 = 0,978. Vagyis a $2$ -nél nagyobb számok biztosan jók, az ennél sokkal kisebbek nem. A pontos határ, amit az octave-ban a norminv(0.975) paranccsal számíthatunk ki (az átrendezés pedig azért érvényes, mert $\Phi$ monoton növő függvény):

$\begin{equation*}\begin{aligned}\frac{0{,}01\sqrt n}{\sqrt{p(1-p)}}&\geq \Phi^{-... ...)=1{,}96\\ n&\geq \frac{1{,}96^2 p(1-p)}{0{,}01^2}. \end{aligned}\end{equation*}$

(8)

Tehát, ha ez a feltétel teljesül, és a közelítések során nem vétettünk túl nagy hibát, akkor ez megadja, hogy legalább hány teszt elvégzésre van szükség. Azonban $p$ -t, ami a védettek valódi aránya, nem ismerjük, pont ezt szeretnénk megbecsülni, így a feltételnek minden $p$ -re teljesülnie kell. Szerencsére $p(1-p)$ nem lehet bármilyen nagy (ez $p$ függvényeként egy parabola, aminek van maximuma), például a számtani-mértani közepek közti egyenlőtlenség alapján $p(1-p)\leq ((p+1-p)/2)^2=1/4$ , és pontos egyenlőség is van, $p=1/2$ esetén. Vagyis a szükséges mintaelemszám, amivel a feltételünk már minden $p$ -re teljesül:

$\displaystyle n\geq\frac{1{,}96^2\cdot 0{,}25}{0{,}01^2}={\bf 9604}.$

3. ábra. Annak valószínűsége, hogy legfeljebb $1\%$ -ot téved a becslés, a valós $p$ védettségi arány függvényében, $n=9604$ teszt esetén

Ezt adja tehát a módszerünk a fent megfogalmazott kérdésre. A közelítés után már csak ekvivalens átalakításokat végeztünk, ezért ha a kiindulás helyes volt, akkor ennél kisebb $n$ nem is lehet jó. Azt is láthatjuk, hogy ha $p$ -ről valamilyen előzetes információnk van, például biztosan tudjuk, hogy legfeljebb $15\%$ , akkor kisebb minta is elég, hiszen ekkor $p(1-p)$ -re az $1/4$ helyett jobb felső becslést írhatunk a (8) egyenlőtlenségben.

Viszont a közelítés miatt mindenképpen szükséges valamiféle ellenőrzés. A 3. ábrán az látható, hogy ha a tesztek száma $n=9604$ , akkor a különböző $p$ értékekre mennyi $\mathbb{P}(\vert X/n-p\vert\leq 0{,}01)$ — vagyis mennyi annak valószínűsége, hogy legfeljebb $1\%$ -ot tévedünk. A binomiális eloszlásból számolt értékek mellett a normális eloszlásből a (7) közelítő képlet eredményét is ábrázoltuk, jól látható, hogy a két görbe egymásra illeszkedik (a binomiális eloszlásnál nem a korábban látott közvetlen számolást, hanem a szoftverek beépített függvényeit alkalmaztuk, ezek a kerekítésekből adódó hibákat a lehető legjobban elkerülik). Az ábra alapján megállapíthatjuk, hogy $n=9604$ esetén valóban tetszőleges $p$ -re legalább $95\%$ valószínűséggel az előírt hibahatáron belüli a becslésünk, viszont ebben az esetben $p=1/2$ esetén egyenlőség lenne, ami azt jelenti, hogy kisebb mintaelemszám már nem elegendő ( $n=9603$ teszt és $p=0{,}5$ esetén annak valószínűsége, hogy legfeljebb $1\%$ -ot tévedünk, „csak” $94{,}9926\%$ ).

6. A relatív hiba, kis p esete és a pontosságtól való függés

Az ábrán azt is érdemes megfigyelni, hogy kicsi vagy egyhez közeli $p$ esetén (ahol pedig egyébként a (4) egyenlet közelítése elromlik), jobbnak tűnik a helyzet, mint amikor $p$ közel van $1/2$ -hez (a szimmetria természetes, hiszen mindegy, hogy a védettek vagy a nem védettek arányára vagyunk kíváncsiak). Azonban vegyük észre, hogy az abszolút eltérésre fogalmaztuk meg a feltételt. Vagyis ebbe az is belefér, hogy $p=0{,}1\%$ helyett a becslésünk $1\%$ legyen – ez legfeljebb $1\%$ -nyi különbség, mégis tízszeres relatív hiba az igazi $p$ -hez képest. Vagyis arra is következtethetünk, hogy kisebb (vagy éppen egyhez közelebbi) $p$ értékeket csak nagyobb mintaelemszámmal tudunk megbecsülni. Ha a $p$ -hez számított relatív hibát adjuk meg, akkor $h=0{,}01$ helyére $r\cdot p$ kerül, és ha visszatérünk a (8) egyenlethez, akkor $1{,}96^2(1-p)/(r^2 p)$ lesz az alsó becslés $n$ -re. Vagyis a szükséges mintaelemszám a $p$ -vel nagyjából fordítottan arányosan tart végtelenhez. Ezt az alábbi táblázat mutatja néhány konkrét érték esetén (itt tehát mindig csak az látszik, hogy ha az adott $p$ -re teljesíteni akarjuk a feltételt, akkor mekkora $n$ -re van szükség):

relatív hiba ( $r$ )	$5\%$	$20\%$	$50\%$
$p=0{,}1$	$13830\ (+1)$	$864\ (-3)$	$139\ (+5)$
$p=0{,}01$	$152128\ (-62)$	$9508 \ (-9)$	$1522\ (-12)$
$p=0{,}001$	$1535103 \ (-136)$	$95944\ (+90)$	$15351\ (+17)$

A táblázatból jól látszik, hogy a relatív hibát alacsonyan tartani kicsi vagy egyhez közeli $p$ esetén jóval nehezebb lehet, mint $1/2$ -hez közelebbi $p$ esetén. A táblázat értékeit a fent bemutatott normális közelítéssel számoltuk, de feltüntettük zárójelben, hogy ez mennyivel tér el a valódi értéktől. Ezzel tehát a közelítésünk pontosságáról is képet kaphatunk. A táblázat és az előző képlet pedig együtt azt is sugallják, hogy ha előzetesen nem tudunk egy biztos alsó becslést $p$ -re (például hogy $p\geq 0{,}01\%$ ), akkor nem is tudunk megfelelő mintaelemszámot mondani a relatív hiba adott korlát alatt tartásához. Ez persze onnan is látható, hogy minél kisebb a $p$ , annál többet kell várnunk arra, hogy akár csak egyetlen védett embert is találjunk. Amíg pedig ilyen nincs, addig csak felső becslést tudunk adni a védettség valószínűségére, alsó becslésre és pontos értékre nem számíthatunk.

4. ábra. A szükséges tesztek száma normális közelítéssel néhány megbízhatósági szint mellett a megengedett abszolút hiba függvényében

A kapott összefüggést általánosan felírva megvizsgálhatjuk, hogy hogyan függ a szükséges mintaelemszám értéke az előírt feltételektől. Az előző számolással teljesen hasonló módon az adódik, hogy ha azt szeretnénk, legalább $q$ valószínűséggel legfeljebb $h$ -t tévedjünk, akkor a fenti (6) közelítéssel a szükséges tesztek száma (visszatérve a $h$ abszolút hibára):

$\displaystyle n\geq \frac{\Phi^{-1}\big(\frac{q+1}{2})}{4\cdot h^2},$

ahol tehát $\Phi(t)=\int_{-\infty}^t \frac{1}{\sqrt{2\pi}}e^{-x^2/2}\,dx$ , és ennek az inverzére van szükség (az octave-ban norminv). Az ebből könnyen látszik, hogy mivel az előírt $h$ hibahatárnak a négyzete szerepel a nevezőben, kétszer akkora pontossághoz négyszer annyi mérés kell. A 4. ábra ezt az összefüggést mutatja néhány különböző $q$ érték esetén. Az ábrán az is látszik, hogy a megbízhatósági szinttől, $q$ -tól való függés kevésbé erős, $99\%$ -os megbízhatósághoz nem lesz szükség nagyságrendekkel több tesztre, mint $95\%$ -oshoz. Persze, ahogy az ábra alapján is sejthető és a képletből következik, az igaz, hogy ha $q\rightarrow 1$ , akkor a szükséges tesztek száma végtelenhez tart: $\Phi^{-1}((q+1)/2)$ az a szám, amitől balra a 2. ábrán $(q+1)/2$ terület van. De mivel a teljes görbe alatti terület $1$ , ez az érték végtelenhez tart, ahogy $q\rightarrow 1$ . Ezek a jelenségek a cikkünk második részében is előjönnek majd, amikor konfidenciaintervallumot adunk meg, vagyis egy olyan intervallumot, ami nagy valószínűséggel tartalmazza a valódi $p$ -t.

7. Elmélet és gyakorlat

A fentiekben tehát a binomiális eloszlás normális közelítésének segítségével meghatároztuk, hogy ahhoz, hogy adott valószínűséggel adott hibahatáron belül legyen a becslésünk, hány teszt szükséges. Végül azt vizsgáljuk meg, hogy a feltevéseink közül mi az, ami a valóságban is könnyen teljesíthető, és mi az, ami további vizsgálatot igényelne.

$\bullet$ Függetlenség. Az egyik feltételezésünk az volt, hogy a tesztek eredménye egymástól független. A vizsgálatba bevont személyek kiválasztásánál ezt könnyű majdnem tökéletesen teljesíteni: ha sorban haladva mindig elfelejtenénk az eddig kiválasztottakat, és úgy választanánk a következő embert, akkor ez teljesülne, ha pedig figyelünk arra, hogy senkit ne válasszunk kétszer, akkor is ez a legtöbb esetben közelítőleg teljesül. Ugyanis mondjuk 2000000 ember közül 10000-t az még könnyen előfordulhat, hogy néhányan többször is szerepelnek, de arra már kevés az esély, hogy valakit olyan sokszor választunk, hogy ettől a mintavételezés eredménye is lényegesen megváltozik a csupa különböző ember esetéhez képest (ez a binomiális eloszlás hipergeometrikus eloszlással való közelítése, [9]). A mérések során ez a feltétel azt jelenti, hogy az egymás után, esetleg párhuzamosan, csoportosítva végzett eljárások eredménye között se legyen összefüggés.

$\bullet$ Azonos védettségi valószínűség. Feltételeztük, hogy minden ember azonos $p$ valószínűséggel védett. Ez egyrészt azt jelenti, hogy a védettek aránya a lakosságban ne változzon jelentősen a tesztek elvégzésének időszaka alatt, ami akár egy-két hét is lehet. Másrészt, ha további tulajdonságokat is figyelembe veszünk (például a lakóhelyet, vagy hogy volt-e az illető mostanában beteg), akkor ez természetesen nem teljesül, de ez nem baj, amíg megelégszünk azzal, hogy az átlagos védettségi valószínűséget vizsgáljuk, és nem az egyes csoportokét külön-külön. Annak valószínűsége ugyanis, hogy egy véletlenszerűen választott ember védett, éppen az átlagos védettségi valószínűség.

$\bullet$ Egyenletes választás. Itt azonban nagyon nem mindegy, hogy a „véletlenszerűen” és az „átlagos” mire vonatkozik. Ideális esetben az átlagban mondjuk egy ország minden lakosa azonos súllyal szerepel. Azt azonban, hogy úgy végezzük a vizsgálatot, hogy mindenkit azonos valószínűséggel válasszunk, legalább annyira nehéz megvalósítani, mint amilyen könnyű leírni. Ha az, hogy valaki védett-e, nem függne attól, hogy milyen könnyű őt egy ilyen vizsálatba bevonni, ez sem jelentene problémát – de leggyakrabban ezt nem tételezhetjük fel. Ennek az úgynevezett mintavételi torzításnak a kiküszöbölésére egy szokásos eljárás az, hogy ha mondjuk a mintában több budapesti ember szerepel, mint amennyi a budapestiek valós aránya, akkor a budapestiek és a vidékiek esetében is külön kiszámítják az átlagot, majd pontosabb (például népszámlálási) adatok alapján ezeket olyan súlyokkal veszik figyelembe, ami ezeknek a csoportoknak a valódi aránya lenne (bővebben például: [5]).

$\bullet$ A tesztek megbízhatósága. Végig úgy számoltunk, hogy a tesztek sosem hibáznak. Ez természetesen a valóságban nincs így, és, főleg kis védettségi valószínűség esetén, még az is nagyon könnyen előfordulhat, hogy a hamis pozitív tesztek száma több, mint valódi pozitívaké. Hiszen ha mondjuk $n=10000$ ember közül csak $50$ védett valójában, akkor a tesztnek körülbelül $1\%$ -os valószínűségű tévedése is elég ahhoz, hogy legyen majdnem $100$ ember, akit hamisan védettnek gondolunk, és hogy a védettség valószínűségét háromszorosan túlbecsüljük. Ez a bizonytalanság, amit információs vagy mérési torzításnak neveznek, mint minden információhiány, azonban feltehetően csak növelné a szükséges tesztek számát, így ebből a szempontból a fenti értékeket alsó becslésnek tekinthetjük. Ennek a hibának a kiküszöbölése az epidemiológiának fontos témaköre, további hibalehetőségek figyelembevételével együtt [10].

$\bullet$ A közelítésből adódó hiba. A fenti levezetésben egy közelítést alkalmaztunk, és nem bizonyítottuk, csak számítógéppel ellenőriztük, hogy a megadott szám valóban elég sok tesztet jelent. Ezen például a már említett Berry–Esséen-tétel vagy ennél is pontosabb becslések alapján lehetne javítani. Ugyanakkor ezekből állításokból is látszik, hogy a normális eloszlással való közelítés elromlik, ha $p$ közel van nullához vagy egyhez. Ilyenkor egy másik lehetőség a Poisson-eloszlással való közelítés lenne [3], [9], de az érvényes marad, hogy ha a relatív hibát is alacsonyan akarjuk tartani, akkor jóval több tesztre lehet szükség. Arra bizonyítást adni, hogy egy adott számú minta már elég a megadott pontossághoz, egyszerűbb módszerekkel, például a Csebisev-egyenlőtlenséggel [2], [3], [9] is lehetséges. Ez azt állítja, hogy annak valószínűsége, hogy egy valószínűségi változó legalább $t$ -vel tér el a saját várható értékétől, legfeljebb a szórásnégyzetének és $t^2$ -nek a hányadosa. Ez a megközelítés azonban olyan egyenlőtlenséget használ, ami minden eloszlásra egyaránt érvényes, és ha speciálisan a binomiális eloszlásra alkalmazzuk, akkor nem túl erős. Konkrétan a $h=1\%$ hiba és $q=95\%$ megbízhatóság mellett ebből az jönne ki, hogy $n\geq 50000$ elég, ami igaz, de azt láttuk, hogy valójában az ötödénél is kevesebb, $9604$ teszt már elég ugyanehhez.

8. Következik: konfidenciaintervallum és maximumlikelihood-becslés

Összességében tehát a fenti matematikai módszerek alkalmazása során egy jól megbízható teszt és a vizsgálatba bevont személyek megfelelően elvégzett kiválasztása mellett számoltunk, közelítő módszerekkel, ahol azonban a közelítés jól működik, amíg a védettek valódi aránya sem nem túl kicsi, sem nem túl nagy a mintaelemszámhoz képest. Ez tehát arra mindenképpen hasznos, hogy lássuk, hogy mi az a szám, aminél semmiképpen nem mehetünk lejjebb az elvégzett tesztek számával adott pontosság mellett. Továbbá azt is megállapíthattuk, hogy a pontosság növelése milyen ütemben növeli a szükséges tesztek számát, és arra is több szempontot láttunk, hogy miért van szükség több tesztre, ha a védettek valódi aránya kicsi, és a relatív hibát akarjuk alacsonyan tartani.

Cikkünk második részében egyrészt arra térünk majd ki, hogy a fenti módszer arra is használható, hogy adott mintaelemszám esetén megfogalmazzuk, hogy mennyire vagyunk biztosak a $p$ -re vonatkozó becslésünkben, amire egy szokásos módszer a konfidenciaintervallum megadása. Ez egy olyan intervallum, ami az ismeretlen $p$ értéket nagy valószínűséggel tartalmazza (bár, mint látni fogjuk, ez nem pontosan azt jelenti, amit elsőre gondolnánk, ugyanis az eddig használt statisztikai felépítést követve $p$ -t továbbra sem tekintjük véletlennek). A fenti módszeren alapuló számításokkal hasonló jelenségeket figyelhetünk majd meg, például hogy megfelelő feltételek mellett fele olyan hosszú konfidenciaintervallumhoz négyszer akkora minta szükséges, vagy hogy kis $p$ esetén a konfidenciaintervallum végpontjainak aránya nagy lehet.

Másrészt egy olyan statisztikai módszert is bemutatunk, amely amellett érvel, hogy minden vizsgált személyhez ugyanazt a $p$ védettségi valószínűséget feltételezve miért az $X/n$ relatív gyakoriság a legjobb becslés, miért nem például $(X+5)/(n+10)$ (amit egy másik módszerrel, a Bayes-becsléssel kaphatnánk, megfelelő előzetes feltételezések mellett). Ezt maximumlikelihood-becslésnek hívják, és a tesztek nyelvén megfogalmazva az alábbi elvből indul ki. Tegyük fel, hogy a laboratóriumban két dobozban vannak tesztek, és arra is emlékszünk, hogy az egyik dobozban a jó tesztek voltak, amik a védettséget $90\%$ valószínűséggel kimutatják, a másikban a rosszak, amik csak $30\%$ valószínűséggel jelzik a védettséget. Csak éppen azt felejtettük el, hogy melyik doboz melyik. Ezért szerzünk néhány biztosan pozitív (védett embertől származó) mintát, és az egyik doboz tíz tesztjével megvizsgáljuk. A tíz mérésből 8 pozitív lett. Biztosak nem lehetünk, de ebben az esetben könnyű megmondani, hogy mire gondolnánk inkább, melyik doboz volt az, amit kipróbáltunk. Ez a gondolat vezet el tehát egy jóval általánosabb statisztikai módszerhez, melyet a következő részben a tesztelés példáján keresztül mutatunk majd be.

Irodalomjegyzék

: [1] Douglas Altman, David Machin, Trevor Bryant, Martin Gardner, Statistics with Confidence: Confidence Intervals and Statistical Guidelines. Second Edition, John Wiley & Sons, New York, 2000.; [2] Arató Miklós, Prokaj Vilmos, Zempléni András, Bevezetés a valószínűségszámításba és alkalmazásaiba: példákkal, szimulációkkal, 2013. https://ttk.elte.hu/dstore/document/901/zempleni.pdf; [3] Csiszár Villő, Valószínűségszámítás 1. http://csvillo.web.elte.hu/mtval/jegyzet.pdf; [4] William Feller, An introduction to probability theory and its applications. Vol. I. Third edition, John Wiley & Sons, Inc., New York, 1968.; [5] David Freedman, Robert Pisani, Roger Purves, Statistics. Fourth edition, W. W. Norton & Company, New York, 2007.; [6] Pierre-Simon de Laplace, Théorie analytiques de probabilités. 3ème éd., revue et augmentée par l'auteur, Courcier, Paris, 1820.; [7] Abraham de Moivre, Miscellanea analytica de seriebus et quadraturis. Lib. 5. J. Tonson & J. Watts, London, 1730.; [8] Rényi Alfréd, Valószínűségszámítás. Tankönyvkiadó, Budapest, 1954.; [9] Sheldon Ross, A first course in probability. Second edition, Macmillan Co., New York, 1984.; [10] Kenneth J. Rothman, Epidemiology: An introduction. Se Oxford University Press, New York, 2012.