Skriv =NORM.DIST(115,100,15,TRUE) inn i en celle og Excel returnerer 0.8413447 uten seremonier. Kallet leses som et oppslag. Det er det ikke. Bak det ene tallet ligger den kumulative normalfordelingen, et integral uten lukket form, og bak CHISQ.INV.RT og BETA.DIST sitter spesielle funksjoner som et nøyaktig bibliotek må evaluere, ikke tilnærme for hånd. En regnearkkomponent som hevder Excel-kompatibilitet må reprodusere disse verdiene til det siste sifferet Excel viser, noe som betyr å reprodusere de numeriske metodene, ikke bare funksjonsnavnene.
HotXLS implementerer mer enn femti av disse statistiske funksjonene, og arbeidet som gjør dem korrekte er nesten helt usynlig fra formellinjen. Dette er en omvisning i hvordan motoren beregner dem: den delte kjernefunksjonen for spesielle funksjoner, forgreningene som holder aritmetikken stabil, og én invers-normal feil som gjemte seg i halen i lang tid fordi det vanlige tilfellet aldri berørte den ødelagte linjen.
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;Metoden Calculate på arbeidsboken kompilerer og evaluerer en ad-hoc-formel mot det aktive arket og returnerer en Variant. En detalj snubler folk i på første forsøk: formelparseren bak Calculate tar semikolon som argumentseparator, så det blir =SUM(A1;B1), ikke =SUM(A1,B1). Lagrede celleformler beholder Excel-standard komma. Den samme evaluatoren sender ut hver statistiske funksjon under, så når en av disse fungerer i Calculate, følger resten samme sti.
De to funksjonene alt annet er bygget på
De fleste kumulative fordelinger i dette settet beregnes ikke ved å summere eller integrere sine egne definisjoner. De beregnes fra to spesielle funksjoner: den regulariserte nedre ufullstendige gammaen, skrevet P(a, x), og den regulariserte ufullstendige betaen, skrevet Ix(a, b). Internt er disse hjelperne som dispatcherne lener seg på, og kjeden er kort. Chi-kvadrat-CDF-en er gamma-CDF-en med form df/2 og skala 2. Gamma-CDF-en er P(a, x) direkte. De kumulative t-, F- og binomialfunksjonene er alle verdier av den regulariserte ufullstendige betaen ved de rette argumentene. Poissons CDF er den øvre ufullstendige gammaen Q. Implementer gamma- og betafunksjonene godt, så arver et dusin fordelinger nøyaktigheten deres gratis.
Ordet "regularisert" er hele poenget. Den rå ufullstendige gammaen vokser som et fakultet, og det rå beta-integralet kan underløpe eller overløpe lenge før svaret gjør det. De regulariserte formene are delt på den fullstendige gammaen eller betaen, så de lever helt i intervallet fra null til én, som er akkurat det området en sannsynlighet opptar. Den normaliseringen er det som lar den samme rutinen betjene en chi-kvadrat med to frihetsgrader og en med to hundre uten at de mellomliggende leddene løper ut av enden på en double. Det forklarer også hvorfor du ikke beregner en CDF bygger på å summere opp en lang hale av tetthetsledd: hvert ledd bærer sin egen avrundingsfeil, feilene akkumuleres etter hvert som serien kjører, og den regulariserte spesialfunksjonen unngår summen komplett ved å evaluere en raskt konvergerende serie eller et kjedebrøk i stedet.
Serier under diagonalen, kjedebrøk over den
Den ufullstendige gammarutinen tar én avgjørelse før den beregner noe: den sammenligner x med a + 1. Den grensen er ikke tilfeldig. Potensserie-ekspansjonen til P(a, x) konvergerer raskt når xer liten i forhold til a, og sakte, til slutt unyttig, når x er stor. Kjedebrøket har motsatt karakter. Så motoren bruker potensserien for x under a + 1 og et Lentz-kjedebrøk for x på eller over a + 1, og hver gren blir bare bedt om å gjøre arbeidet den er god på.
Kjedebrøket trenger én beskyttelse. Lentz' metode fungerer ved å bære en løpende teller og nevner og invertere nevneren på hvert trinn, og hvis en av dem nærmer seg null, sprekker inverteringen. Løsningen er en bitteliten bunngrense: når et mellomledd faller under omtrent 1e-30 i størrelse, låses det til 1e-30, noe som holder rekursjonen endelig uten å forstyrre den konvergerte verdien. Den samme låsen vises i den ufullstendige betaens kjedebrøk av samme grunn. Det er en liten konstant som gjør viktig arbeid, forskjellen mellom en stabil evaluering og en divisjon med noe som ikke kan skyves fra null.
Den øvre halen, Q(a, x), er ganske enkelt 1 minus P(a, x), og det er slik Poissons kumulative gren beregnes: sannsynligheten for maksimalt k hendelser med gjennomsnitt λ er Q(k + 1, λ). Routing it through the upper incomplete gamma rather than summing k + 1 Poisson terms is, again, a choice to evaluate one convergent expression instead of accumulating many small ones.
Diskrete masser uten fakultetsoverløp
De diskrete fordelingene reiser en annen fare. En binomial sannsynlighetsmasse involverer en binomialkoeffisient, og koeffisienten for "52 velg 26" er et enormt heltall. Oppretter du det direkte, vil telleren overløpe en double før divisjonen som ville bringe det tilbake til en fornuftig sannsynlighet. Motoren oppretter det aldri direkte. Den beregner fakultetene i log-rommet via log-gamma-funksjonen, adderer og subtraherer loggene, bretter inn loggen for suksess- og fiaskosannsynlighetene, og tar eksponenten én gang helt til slutt.
// 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));Selve log-gamma-funksjonen er en Lanczos-tilnærming, nøyaktig over hele den positive aksen og billig å evaluere. Fordi hver store megde holdes som sin logaritme frem til den endelige Exp, det største tallet rutinen noensinne materialiserer selve sannsynligheten, som er maksimalt én. Poissons massefunksjon følger samme oppskrift, der det enkle log-gamma-leddet står i stedet for fakultetet i nevneren. De lukkede formene er spesialbehandlet i kantene, der p er nøyaktig null eller én, slik at koden anvender aldri kaller Ln(0). HotXLS returnerer 0.2460938 for BINOM.DIST(5,10,0.5,FALSE) og 0.6766764 for den kumulative POISSON.DIST(2,2,TRUE), noe som samsvarer med Excel gjennom sifrene den skriver ut.
Inverser ved å avgrense den fremadgående kurven
En invers fordelingsfunksjon stiller det motsatte spørsmålet: gitt en sannsynlighet, finn verdien der CDF-en er lik den. Bare én invers i dette settet har a rask direkte formel. NORM.S.INV, den inverse standardnormalen, bruker en Acklam-rasjonell tilnærming, et par polynomske forhold nøyaktig til omtrent presisjonen av en double over hele området, delt inn i en sentral region og to haler. Det er en evaluering på lukket form uten iterasjon.
De andre inversene har ingen slik formel, so motoren inverterer dem numerisk. Den avgrenser (brackets) svaret med en nedre og øvre grense valgt fra fordelingens støtte, og halverer (bisects) deretter: evaluer den fremadgående CDF-en ved midtpunktet, flytt den grensen som holder målsannsynligheten innesluttet, og gjenta til intervallet er smalt. For gamma- og chi-kvadrat-inversene starter grensen på null og et generøst øvre estimat bygget ut fra formen og skalaen, og dobler den øvre grensen hvis sannsynligheten ennå ikke er innesluttet. Den inverse t-en avgrenser symmetriske grenser som utvides utover; den inverse F-en halverer på et ikke-negativt intervall. Kostnaden er noen titalls CDF-evalueringer per kall, noe som er usynlig ved regnearkhastighet, og fordelen er at hver invers er nøyaktig like nøyaktig som den fremadgående funksjonen den inverterer. Det er derfor en rundreise som CHISQ.DIST(CHISQ.INV(0.7,5),5,TRUE) returnerer 0.7 innenfor et hårstrå.
Base-10-logaritmen som gjemte seg i halen
Her er feilen som er verdt å fortelle, fordi det er typen som overlever lenge. Den inverse normalrutinen til Acklam har tre grener. Den brede sentrale grenen, som brukes når sannsynligheten ligger mellom omtrent 0.025 og 0.975, kjører inndataene gjennom et polynomforhold uten logaritme noen steder i seg. De to halegrenene, for svært små eller svært store sannsynligheter, tar hver en logaritme av inndataene først, fordi halen oppfører seg som kvadratroten av minus den naturlige logaritmen til p.
En tidlig versjon av halegrenen tok en 10-er-logaritme (base-10) der den naturlige logaritmen hørte hjemme. De to skiller seg med en konstant faktor på omtrent 2.30, så haleresultatene var feil med en konsistent, betydelig margin. Og yet funksjonen så fin ut i enhver tilfeldig sjekk, fordi tilfeldige sjekker lever i midten. NORM.S.INV(0.5) er null, NORM.S.INV(0.975) er læreboken 1.959964, og begge disse kjører gjennom det sentrale polynomet som aldri kaller en logaritme i det hele tatt. Feilen dukket bare opp når en sannsynlighet krysset inn i en hale, for eksempel NORM.S.INV(0.001), som må returnere -3.0902323, men i stedet kom tilbake skjevt på grunn av forholdet mellom naturlig logaritme og 10-er-logaritme. Enhver funksjon som avhenger av den inverse normalen i halen sin, inkludert konfidensintervall-hjelperne, arvet den samme skjevheten. Leksjonen er hverdagslig og kostbar: en funksjon med en grenstruktur trenger testpunkter inne i hver gren, fordi en korrekt vanlig sti gladelig vil maskere en ødelagt sjelden en. Løsningen var en endring av ett token fra base-10 logg til naturlig logg, og haleverdiene spratt på plass til Excels.
Tegnet på x bestemmer t-fordelingens hale
Student t-kumulativfunksjonen har en finurlighet som er lett å ta baklengs. Verdien kommer fra den regulariserte ufullstendige betaen evaluert ved df / (df + x²), men den betaverdien er sannsynligheten i halen utenfor størrelsen på x, not den kumulative sannsynligheten opp til x. Den symmetriske formen til t-fordelingen betyr at konverteringen avhenger av hvilken side av null x faller på.
// 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 meanIngenting av dette er synlig fra cellen, og det er det tilsiktede resultatet. Du skriver en formel og leser et tall som stemmer overens med Excel. Hvis du utvider motoren med din egen logikk, er mekanikken for å registrere en funksjon dekket i vår gjennomgang av formelmotoren og egendefinerte funksjoner, og måten formler når på tvers av ark og navngitte områder er dekket i artikkelen om definerte navn og formler på tvers av ark. Alt dette leveres inne i HotXLS regnearkkomponent for Delphi og C++Builder, sammen med lesing, skriving, diagrammering og formaterings-API-ene som er dekket andre steder på denne bloggen.