Upišite =NORM.DIST(115,100,15,TRUE) u ćeliju i Excel bez puno objašnjavanja vraća 0.8413447. Taj poziv zvuči kao obično pretraživanje, ali to nije. Iza tog jednog broja krije se kumulativna normalna raspodjela, integral bez zatvorenog oblika, a iza CHISQ.INV.RT i BETA.DIST nalaze se posebne funkcije koje pažljivo izgrađena knjižnica mora točno izračunati, a ne samo ručno procijeniti. Komponenta proračunske tablice koja tvrdi da je kompatibilna s Excelom mora reproducirati ove vrijednosti do zadnje znamenke koju Excel prikazuje, što znači reproducirati numeričke metode, a ne samo nazive funkcija
HotXLS implementira više od pedeset takvih statističkih funkcija, a rad koji ih čini točnima gotovo je potpuno nevidljiv iz trake formule. Ovo je pregled načina na koji ih pokretač izračunava: zajednička jezgra posebnih funkcija, odluke o granama koje aritmetiku drže stabilnom, te jedna pogreška inverzne normalne raspodjele koja se dugo skrivala u repu jer uobičajeni slučajevi nikada nisu dotaknuli taj neispravni kôd
Jedan poziv u radnom listu, a iza njega pedeset raspodjela
Funkcije obuhvaćaju obitelji koje se često koriste u statističkim radnim knjigama. Tu je normalna obitelj, NORM.DIST i NORM.S.DIST s njihovim inverznim funkcijama; obitelj gama i hi-kvadrat, GAMMA.DIST, CHISQ.DIST, CHISQ.DIST.RT, CHISQ.INV.RT; obitelj beta, BETA.DIST i BETA.INV; raspodjele uzoraka T.DIST, T.DIST.2T, F.DIST i F.INV; diskretni par BINOM.DIST i POISSON.DIST; te pomoćnici za zaključivanje kao što su CONFIDENCE.T i CONFIDENCE.NORM. Iz perspektive pozivatelja, svaka od njih je jednostavna formula. Postavite ulaze u ćelije, zatražite od radne knjige izračun i pročitate rezultat
Metoda Calculate na radnoj knjizi prevodi i procjenjuje ad-hoc formulu u odnosu na aktivni list i vraća Variant. Jedan detalj može zbuniti programere pri prvom pokušaju: raščlanjivač formula iza metode Calculate koristi točku-zarez kao separator argumenata, pa se piše =SUM(A1;B1), a ne =SUM(A1,B1). Formule pohranjene u ćelijama zadržavaju standardni zarez iz Excela. Isti procjenitelj prosljeđuje svaku statističku funkciju u nastavku, pa čim jedna od njih proradi u metodi Calculate, ostale slijede istu putanju
Dvije funkcije na kojima se gradi sve ostalo
Većina kumulativnih raspodjela u ovom skupu ne izračunava se zbrajanjem ili integriranjem vlastitih definicija. Izračunavaju se iz dviju posebnih funkcija: regularizirane donje nepotpune gama funkcije, zapisane kao P(a, x), i regularizirane nepotpune beta funkcije, zapisane kao Ix(a, b). Interno, to su pomoćne funkcije na koje se prosljeđivači oslanjaju, a lanac je kratak. CDF (kumulativna funkcija raspodjele) hi-kvadrata je gama CDF s oblikom df/2 i skalom 2. Gama CDF je izravno P(a, x). Kumulativne funkcije za t, F i binomnu raspodjelu predstavljaju vrijednosti regularizirane nepotpune bete s odgovarajućim argumentima. Poissonov CDF je gornja nepotpuna gama Q. Dobro implementirajte gama i beta funkcije i desetak drugih raspodjela naslijedit će njihovu točnost bez ikakvog dodatnog truda
Riječ "regularizirana" ključ je cijele priče. Sirova nepotpuna gama raste poput faktorijela, a sirovi beta integral može doživjeti podlijevanje (underflow) ili prelijevanje (overflow) puno prije samog odgovora. Regularizirani oblici podijeljeni su potpunom gamom ili betom, pa u potpunosti žive u intervalu od nule do jedan, što je točno raspon koji vjerojatnost zauzima. Ta normalizacija je ono što omogućuje da ista rutina služi za hi-kvadrat s dva stupnja slobode i za onaj s dvije stotine, a da međuvrijednosti ne pobjegnu izvan granica tipa double. To također objašnjava zašto kumulativnu funkciju ne računate zbrajanjem dugog repa gustoće: svaki član nosi vlastitu pogrešku zaokruživanja, pogreške se nakupljaju kako se niz nastavlja, a regularizirana posebna funkcija potpuno zaobilazi zbroj procjenom brzo konvergirajućeg niza ili verižnog razlomka
Niz ispod dijagonale, verižni razlomak iznad nje
Rutina za nepotpunu gamu donosi jednu odluku prije nego što bilo što izračuna: uspoređuje x s a + 1. Ta granica nije proizvoljna. Razvoj u potencijalni red za P(a, x) brzo konvergira kada je x malen u odnosu na a, a sporo, na kraju i beskorisno, kada je x velik. Verižni razlomak ima suprotne karakteristike. Stoga pokretač koristi potencijalni red za x manji od a + 1 i Lentzov verižni razlomak za x jednak ili veći od a + 1, tako da se u svakoj grani obavlja samo onaj posao za koji je najprikladnija
Verižni razlomak zahtijeva jednu zaštitu. Lentzova metoda radi tako da prenosi tekući brojnik i nazivnik i invertira nazivnik pri svakom koraku, a ako se bilo koji od njih približi nuli, inverzija ne uspijeva. Rješenje je vrlo niska donja granica: kad god magnitude nekog međuvremena padne ispod otprilike 1e-30, ona se postavlja (clamp) na 1e-30, što održava rekurziju konačnom bez narušavanja konvergirane vrijednosti. Ista se granica pojavljuje i u verižnom razlomku nepotpune bete iz istog razloga. To je mala konstanta koja obavlja važan zadatak — razlika između stabilne procjene i dijeljenja s nečim što se ne razlikuje od nule
Gornji rep, Q(a, x), jednostavno je 1 minus P(a, x), i tako se računa Poissonova kumulativna grana: vjerojatnost od najviše k događaja sa srednjom vrijednošću λ je Q(k + 1, λ). Usmjeravanje kroz gornju nepotpunu gamu umjesto zbrajanja k + 1 Poissonovih članova je, ponovno, odluka da se procijeni jedan konvergirajući izraz umjesto akumuliranja mnogih malih
Diskretne mase bez prelijevanja faktorijela
Diskretne raspodjele donose drugačiju opasnost. Binomna masa vjerojatnosti uključuje binomni koeficijent, a koeficijent za "pedeset i dva povrh dvadeset i šest" je ogroman cijeli broj. Formirajte ga izravno i brojnik će prelijevati tip double prije dijeljenja koje bi ga vratilo na razumnu vjerojatnost. Pokretač ga nikada ne stvara izravno. On izračunava faktorijele u logaritamskom prostoru putem log-gama funkcije, zbraja i oduzima logaritme, dodaje logaritme vjerojatnosti uspjeha i neuspjeha, te provodi potenciranje (eksponencijalni račun) tek na samom kraju
// Binomial probability mass, evaluated entirely in log space.
// LnGammaF(n+1) is ln(n!); the three log-factorials form ln(C(n,k)),
// and the whole exponent is built before a single Exp call.
// ln P(X=k) = ln(n!) - ln(k!) - ln((n-k)!) + k*ln(p) + (n-k)*ln(1-p)
result := Exp(LnGammaF(nt + 1) - LnGammaF(kk + 1) - LnGammaF(nt - kk + 1)
+ kk * Ln(pp) + (nt - kk) * Ln(1 - pp));Sama log-gama funkcija je Lanczoseva aproksimacija, točna na cijeloj pozitivnoj osi i jeftina za procjenu. Budući da se svaka velika količina drži kao svoj logaritam do konačnog Exp, najveći broj koji rutina ikada stvori je sama vjerojatnost, koja iznosi najviše jedan. Poissonova funkcija mase slijedi isti recept, s jednim log-gama članom koji zamjenjuje faktorijel u nazivniku. Zatvoreni oblici se posebno rješavaju na rubovima, gdje je p točno nula ili jedan, pa kôd nikada ne poziva Ln(0). HotXLS vraća 0.2460938 za BINOM.DIST(5,10,0.5,FALSE) i 0.6766764 za kumulativnu POISSON.DIST(2,2,TRUE), što odgovara Excelu u svim prikazanim znamenkama
Inverzne vrijednosti uokvirivanjem izravne krivulje
Inverzna funkcija raspodjele postavlja suprotno pitanje: za zadanu vjerojatnost, pronađite vrijednost čiji joj je CDF negativan. Samo jedna inverzna funkcija u ovom skupu ima brzu izravnu formulu. NORM.S.INV, inverzna standardna normalna raspodjela, koristi Acklamovu racionalnu aproksimaciju, par omjera polinoma točnih na razini preciznosti tipa double kroz cijeli raspon, podijeljen na središnju regiju i dva repa. To je procjena zatvorenog oblika bez iteracije
Ostale inverzne funkcije nemaju takvu formulu, pa ih pokretač invertira numerički. On uokviruje (brackets) odgovor s niskom i visokom granicom odabranom iz podrške raspodjele, a zatim dijeli na pola (bisekcija): procjenjuje izravni CDF na središnjoj točki, pomiče onu granicu koja zadržava ciljnu vjerojatnost zatvorenom i ponavlja postupak dok interval ne postane uzak. Za inverzne vrijednosti game i hi-kvadrata, okvir počinje od nule i velikodušne gornje procjene izgrađene na temelju oblika i skale, udvostručujući gornju granicu ako vjerojatnost još nije uokvirena. Inverzni t uokviruje simetrične granice koje se šire prema van; inverzni F dijeli na pola na nenegativnom intervalu. Cijena je nekoliko desetaka CDF procjena po pozivu, što je neprimjetno pri brzini rada proračunske tablice, a prednost je što je svaka inverzna funkcija točno onoliko točna koliko i izravna funkcija koju invertira. Zato dvosmjerni proces poput CHISQ.DIST(CHISQ.INV(0.7,5),5,TRUE) vraća 0.7 gotovo u dlaku točno
Baza-10 logaritam koji se skrivao u repu
Evo pogreške koju vrijedi ispričati, jer se radi o vrsti koja dugo opstaje. Acklamova inverzna normalna rutina ima tri grane. Široka središnja grana, koja se koristi kad god vjerojatnost leži između otprilike 0.025 i 0.975, provlači ulaz kroz polinomski omjer bez ikakvog logaritma. Dvije grane na repu, za vrlo male ili vrlo velike vjerojatnosti, najprije uzimaju logaritam ulaza, jer se rep ponaša kao kvadratni korijen minus prirodnog logaritma od p
Rana verzija grane na repu uzimala je logaritam s bazom 10 tamo gdje je pripadao prirodni logaritam. Njih se dvoje razlikuju za konstantan faktor od oko 2.30, pa su rezultati u repu bili pogrešni za dosljednu, znatnu marginu. Pa ipak, funkcija je izgledala sasvim u redu pri svakoj usputnoj provjeri, jer te provjere obično obitavaju u sredini. NORM.S.INV(0.5) je nula, NORM.S.INV(0.975) je udžbenički 1.959964, a oba ta slučaja prolaze kroz središnji polinom koji uopće ne poziva logaritam. Pogreška se pojavljivala tek kada bi vjerojatnost prešla u rep, na primjer NORM.S.INV(0.001), koji mora vratiti -3.0902323, a umjesto toga se vraćao s odstupanjem zbog omjera prirodnog logaritma i logaritma s bazom 10. Svaka funkcija koja ovisi o inverznoj normalnoj raspodjeli u svom repu, uključujući pomoćnike za intervale pouzdanosti, naslijedila je isto odstupanje. Pouka je jednostavna i skupa: funkcija sa strukturom grana treba testne točke unutar svake pojedine grane, jer će ispravna uobičajena putanja veselo maskirati pokvarenu rijetku. Rješenje je bila izmjena jednog simbola iz logaritma s bazom 10 u prirodni logaritam, nakon čega su se vrijednosti u repu odmah uskladile s Excelovima
Predznak od x određuje rep t-raspodjele
Kumulativna funkcija Studentove t-raspodjele nosi suptilnost koju je lako pogrešno shvatiti. Njezina vrijednost dolazi iz regularizirane nepotpune bete procijenjene na df / (df + x²), ali ta beta vrijednost je vjerojatnost u repu izvan magnitude od x, a ne kumulativna vjerojatnost do x. Simetrični oblik t-raspodjele znači da pretvorba ovisi o tome s koje strane nule leži x
// Student t CDF. ib is the regularized incomplete beta at df/(df+x*x),
// which measures the symmetric tail. The cumulative value depends on
// the sign of x; returning ib unconverted gives the wrong tail.
ib := BetaIF(df / 2, 0.5, df / (df + x * x));
if x > 0 then
result := 1 - 0.5 * ib // above the mean: one minus half the tail
else if x < 0 then
result := 0.5 * ib // below the mean: half the tail
else
result := 0.5; // exactly at the meanZa x veći od nule, kumulativna vjerojatnost je jedan minus pola simetričnog repa; za x manji od nule, to je pola tog repa; na nuli je to točno jedna polovina. Vratite li vrijednost bete izravno, prijavit ćete pogrešnu stranu raspodjele, promašivši cijelo tijelo krivulje za bilo koji x različit od nule. Varijante desnog repa i dvaju repova grade se na istoj grani, zbog čega se T.DIST.2T(1,1) vraća kao 0.5, a T.DIST(1,1,TRUE) kao 0.75, dok se inverzni T.INV dijeli na pola u odnosu na ovaj ispravljeni CDF kako bi se dvosmjerni proces uspješno zatvorio
Ništa od ovoga nije vidljivo iz ćelije, i to je željeni ishod. Zapišete formulu i pročitate broj koji se slaže s Excelom. Ako proširujete pokretač vlastitom logikom, mehanika registracije funkcije opisana je u našem vodiču kroz mehanizam formula i prilagođene funkcije, a način na koji formule dopiru do drugih listova i imenovanih raspona opisan je u članku o definiranim imenima i referencama među listovima. Sve to dolazi unutar komponente HotXLS spreadsheet component za Delphi i C++Builder, uz API-je za čitanje, pisanje, grafikone i oblikovanje o kojima se raspravlja drugdje na ovom blogu