Ukucajte =NORM.DIST(115,100,15,TRUE) u ćeliju i Excel vraća 0.8413447 bez ceremonije. Poziv se čita kao obično traženje vrednosti. Ali nije. Iza tog jednog broja stoji kumulativna normalna raspodela, integral koji nema zatvoreni oblik, a iza CHISQ.INV.RT i BETA.DIST leže specijalne funkcije koje pažljiva biblioteka mora da proceni, a ne da aproksimira ručno. Komponenta tabele koja pretenduje na kompatibilnost sa Excel-om mora da reprodukuje ove vrednosti do poslednje cifre koju Excel prikazuje, što znači reprodukovanje numeričkih metoda, a ne samo imena funkcija.
HotXLS implementira više od pedeset ovih statističkih funkcija, a rad koji ih čini ispravnim je skoro potpuno nevidljiv sa trake za formule. Ovo je pregled načina na koji ih mehanizam računa: deljeno jezgro specijalnih funkcija, odluke o grananju koje održavaju aritmetiku stabilnom i jedna greška inverzne normale koja se dugo krila u repu jer uobičajeni slučaj nikada nije dodirnuo pokvarenu liniju koda.
Jedan poziv u radnom listu, pedeset raspodela iza njega
Funkcije obuhvataju porodice koje statistička radna sveska traži. Tu je porodica normalne raspodele, NORM.DIST i NORM.S.DIST sa njihovim inverzama; porodica gama i hi-kvadrat raspodele, GAMMA.DIST, CHISQ.DIST, CHISQ.DIST.RT, CHISQ.INV.RT; porodica beta raspodele, BETA.DIST i BETA.INV; raspodele uzoraka T.DIST, T.DIST.2T, F.DIST i F.INV; diskretni par BINOM.DIST i POISSON.DIST; i pomoćnici za zaključivanje kao što su CONFIDENCE.T i CONFIDENCE.NORM. Sa mesta pozivaoca, svaka je jedna formula. Postavite ulaze u ćelije, zatražite od radne sveske da ih proceni i pročitate rezultat.
var
wb: IXLSWorkbook;
sh: IXLSWorksheet;
begin
wb := TXLSWorkbook.Create;
sh := wb.Sheets.Add;
sh.Range['A1', 'A1'].Value := 115; // observation
sh.Range['A2', 'A2'].Value := 100; // mean
sh.Range['A3', 'A3'].Value := 15; // standard deviation
// The XLS formula parser uses ';' as the argument separator.
Writeln(wb.Calculate('=NORM.DIST(A1;A2;A3;TRUE())')); // 0.8413447
Writeln(wb.Calculate('=CHISQ.INV.RT(0.05;10)')); // 18.3070381
Writeln(wb.Calculate('=BETA.DIST(0.5;2;3;TRUE())')); // 0.6875
end;
Metoda Calculate na radnoj svesci kompajlira i procenjuje ad-hoc formulu u odnosu na aktivni list i vraća Variant. Jedan detalj sapliće korisnike pri prvom pokušaju: parser formula iza metode Calculate uzima tačka-zarez kao separator argumenata, tako da je =SUM(A1;B1), a ne =SUM(A1,B1). Sačuvane formule ćelija zadržavaju Excel-standardni zarez. Isti evaluator otprema svaku statističku funkciju ispod, tako da kada jedna od njih proradi u Calculate, ostale prate istu putanju.
Dve funkcije na kojima je sve ostalo izgrađeno
Većina kumulativnih raspodela u ovom skupu se ne računa sabiranjem ili integracijom sopstvenih definicija. One se računaju iz dve specijalne funkcije: regularizovane donje nepotpune gama funkcije, zapisane kao P(a, x), i regularizovane nepotpune beta funkcije, zapisane kao Ix(a, b). Interno, to su pomoćnici na koje se otpremnici oslanjaju, i lanac je kratak. CDF hi-kvadrata je CDF gama sa oblikom df/2 i skalom 2. Gama CDF je direktno P(a, x). Kumulativne funkcije t, F i binomne raspodele su sve vrednosti regularizovane nepotpune beta funkcije pri odgovarajućim argumentima. Poisson-ov CDF je gornja nepotpuna gama Q. Implementirajte dobro gama i beta funkcije i desetak raspodela će naslediti njihovu tačnost besplatno.
Reč „regularizovana” je cela poenta. Sirova nepotpuna gama raste kao faktorijel, a sirovi beta integral može da pretrpi podlivanje ili prelivanje mnogo pre nego što to učini odgovor. Regularizovani oblici su podeljeni sa potpunom gama ili beta funkcijom, tako da u potpunosti žive u intervalu od nule do jedan, što je tačno opseg koji zauzima verovatnoća. Ta normalizacija je ono što omogućava da ista rutina služi za hi-kvadrat sa dva stepena slobode i za onaj sa dve stotine bez da srednji članovi pobegnu sa kraja opsega dvostruke preciznosti. To takođe objašnjava zašto ne računate CDF sabiranjem dugog repa članova gustine: svaki član nosi sopstvenu grešku zaokruživanja, greške se akumuliraju kako serija teče, a regularizovana specijalna funkcija u potpunosti zaobilazi sumu procenjujući umesto toga brzo konvergentnu seriju ili verižni razlomak.
Serija ispod dijagonale, verižni razlomak iznad nje
Rutina nepotpune gama funkcije donosi jednu odluku pre nego što bilo šta izračuna: poredi x sa a + 1. Ta granica nije proizvoljna. Razvoj u stepeni red funkcije P(a, x) brzo konvergira kada je x malo u odnosu na a, a sporo, na kraju i beskorisno, kada je x veliko. Verižni razlomak ima suprotan karakter. Zato mehanizam koristi stepeni red za x ispod a + 1 i Lentz-ov verižni razlomak za x na ili iznad a + 1, pri čemu se od svake grane traži da obavi samo onaj posao u kom je dobra.
Verižnom razlomku je potrebna jedna zaštita. Lentz-ova metoda radi tako što nosi tekući brojilac i imenilac i invertuje imenilac na svakom koraku, a ako se bilo koji od njih približi nuli, inverzija propada. Rešenje je sićušni donji prag: kad god srednji član padne ispod otprilike 1e-30 u magnitudi, on se fiksira na 1e-30, što održava rekurziju konačnom bez narušavanja konvergirane vrednosti. Isti prag se pojavljuje i u verižnom razlomku nepotpune beta funkcije iz istog razloga. To je mala konstanta koja obavlja noseći rad, čineći razliku između stabilne procene i deljenja nečim što se ne razlikuje od nule.
Gornji rep, Q(a, x), jeste jednostavno 1 minus P(a, x), i tako se računa kumulativna grana Poisson-ove raspodele: verovatnoća od najviše k događaja sa srednjom vrednošću λ jeste Q(k + 1, λ). Njeno usmeravanje kroz gornju nepotpunu gama funkciju umesto sabiranja k + 1 Poisson-ovih članova je, opet, izbor da se proceni jedan konvergentan izraz umesto da se akumulira mnogo malih.
Diskretne mase bez prelivanja faktorijela
Diskretne raspodele donose drugačiju opasnost. Binomna masa verovatnoće uključuje binomni koeficijent, a koeficijent za izbor dvadeset šest od pedeset dva je ogroman ceo broj. Ako ga formirate direktno, brojilac preliva opseg dvostruke preciznosti pre deljenja koje bi ga vratilo na razumnu verovatnoću. Mehanizam ga nikada ne formira direktno. On računa faktorijele u logaritamskom prostoru kroz log-gama funkciju, sabira i oduzima logaritme, uključuje logaritme verovatnoća uspeha i neuspeha i stepenuje samo jednom 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 Lanczos-ova aproksimacija, tačna na celoj pozitivnoj osi i jeftina za procenu. Pošto se svaka velika količina drži kao njen logaritam do konačnog Exp, najveći broj koji rutina ikada materijalizuje je sama verovatnoća, koja iznosi najviše jedan. Poisson-ova funkcija mase prati isti recept, sa jednim log-gama članom koji zamenjuje faktorijel u imeniocu. Zatvoreni oblici se tretiraju kao posebni slučajevi na ivicama, gde je p tačno nula ili jedan, tako da kod nikada ne poziva Ln(0). HotXLS vraća 0.2460938 for BINOM.DIST(5,10,0.5,FALSE) i 0.6766764 za kumulativnu POISSON.DIST(2,2,TRUE), što se podudara sa Excel-om kroz sve cifre koje štampa.
Inverze uokvirivanjem krive unapred
Funkcija inverzne raspodele postavlja suprotno pitanje: za datu verovatnoću, pronaći vrednost čiji je CDF jednak njoj. Samo jedna inverza u ovom skupu ima brzu direktnu formulu. NORM.S.INV, inverzna standardna normalna raspodela, koristi Acklam-ovu racionalnu aproksimaciju, par odnosa polinoma tačnih do otprilike preciznosti dvostruke preciznosti na celom opsegu, podeljenih na centralni region i dva repa. To je procena u zatvorenom obliku bez iteracija.
Ostale inverze nemaju takvu formulu, pa ih mehanizam invertuje numerički. On uokviruje odgovor donjom i gornjom granicom izabranom iz podrške raspodele, a zatim vrši bisekciju: proceni CDF unapred u srednjoj tački, pomeri granicu koja drži ciljnu verovatnoću zatvorenom i ponavlja postupak dok interval ne postane uzak. Za inverze gama i hi-kvadrat raspodele, uokviravanje počinje od nule i velikodušne gornje procene izgrađene na osnovu oblika i skale, duplirajući gornju granicu ako verovatnoća još uvek nije zatvorena. Inverza t-raspodele uokviruje simetrične granice koje se šire ka spolja; inverza F-raspodele vrši bisekciju na ne-negativnom intervalu. Cena je nekoliko desetina CDF procena po pozivu, što je neprimetno pri brzini rada sa tabelama, a korist je to što je svaka inverza tačna upravo onoliko koliko i funkcija unapred koju invertuje. Zato povratno putovanje kao što je CHISQ.DIST(CHISQ.INV(0.7,5),5,TRUE) vraća 0.7 u dlaku tačno.
Logaritam sa osnovom 10 koji se krio u repu
Evo greške koju vredi pomenuti, jer je to vrsta greške koja dugo opstaje. Acklam-ova rutina za inverznu normalnu raspodelu ima tri grane. Široka centralna grana, koja se koristi kad god se verovatnoća nalazi između otprilike 0.025 i 0.975, provlači ulaz kroz polinomijalni odnos bez ikakvog logaritma u sebi. Dve repne grane, za veoma male ili veoma velike verovatnoće, prvo uzimaju logaritam ulaza, jer se rep ponaša kao kvadratni koren iz minus prirodnog logaritma od p.
Rana verzija repne grane je uzimala logaritam sa osnovom 10 tamo gde je pripadao prirodni logaritam. Ova dva se razlikuju za konstantan faktor od oko 2.30, tako da su rezultati u repu bili pogrešni za doslednu, priličnu marginu. Pa ipak, funkcija je izgledala u redu u svakoj usputnoj proveri, jer usputne provere žive u sredini. NORM.S.INV(0.5) je nula, NORM.S.INV(0.975) je školskih 1.959964, a oba ta primera prolaze kroz centralni polinom koji uopšte ne poziva logaritam. Greška se pojavljivala tek kada bi verovatnoća prešla u rep, na primer NORM.S.INV(0.001), koji mora da vrati -3.0902323, a umesto toga se vraćao sa odstupanjem zbog odnosa prirodnog i dekadnog logaritma. Svaka funkcija koja zavisi od inverzne normale u svom repu, uključujući pomoćnike za intervale poverenja, nasledila je isto odstupanje. Ispravka je bila promena jednog tokena sa dekadnog logaritma na prirodni logaritam, i vrednosti u repu su se trenutno uskladile sa Excel-ovim.
Znak od x odlučuje o repu t-raspodele
Kumulativna funkcija Student-ove t-raspodele nosi suptilnost koju je lako shvatiti obrnuto. Njena vrednost dolazi iz regularizovane nepotpune beta funkcije procenjene na df / (df + x²), ali ta beta vrednost je verovatnoća u repu izvan magnitude od x, a ne kumulativna verovatnoća do x. Simetrični oblik t-raspodele znači da konverzija zavisi od toga sa koje strane nule x pada.
// 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 mean
Za x iznad nule, kumulativna verovatnoća je jedan minus polovina simetričnog repa; za x ispod nule, to je polovina tog repa; na nuli, to je tačno jedna polovina. Ako vratite beta vrednost direktno, prijavljujete pogrešnu stranu raspodele, promašivši za celo telo krive za bilo koje ne-nulto x. Varijante desnog repa i sa dva repa se grade 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, a inverzna T.INV vrši bisekciju u odnosu na ovaj ispravljeni CDF tako da se krug zatvara.
Ništa od ovoga nije vidljivo iz ćelije, i to je željeni ishod. Vi napišete formulu i pročitate broj koji se slaže sa Excel-om. Ako proširujete mehanizam sopstvenom logikom, mehanika registrovanja funkcije pokrivena je u našem vodiču kroz mehanizam formula i prilagođene funkcije, a način na koji formule dosežu preko listova i imenovanih opsega pokriven je u članku o definisanim imenima i formulama između listova. Sve to se isporučuje unutar HotXLS komponente za tabele za Delphi i C++Builder, zajedno sa API-jima za čitanje, pisanje, iscrtavanje grafikona i formatiranje koji su pokriveni na drugim mestima na ovom blogu.