Function PoissonInv(Prob As Double, Mean As Double) As Variant
' shg 2011, 2012, 2014, 2015-0415
' For a Poisson process with mean Mean, returns a three-element array:
' o The smallest integer N such that POISSON(N, Mean, True) >= Prob
' o The CDF for N-1 (which is < Prob)
' o The CDF for N (which is >= Prob)
' E.g., POISSON(5, 10, TRUE) returns 0.067085962
' PoissonInv(0.067085962, 10) returns 5
' Returns a descriptive error if Prob <= 0 or Prob >= 1
' If Mean >= 100, then uses a normal approximation
Dim N As Long ' number of events
Dim CDF As Double ' cumulative distribution function at N
Dim CDFOld As Double ' CDF at N-1
' These two variables are used to simplify the probability mass
' function summation CDF = CDF + Exp(-Mean) * Mean ^ N / N!
Dim dExpMean As Double ' =Exp(-Mean)
Dim dK As Double ' incremental power & factorial
If Prob <= 0 Or Prob >= 1 Then
PoissonInv = "Prob ]0,1["
ElseIf Mean < 100 Then
dExpMean = Exp(-Mean)
dK = 1#
CDF = dExpMean
Do While CDF < Prob - 0.000000000000001
CDFOld = CDF
N = N + 1
dK = dK * Mean / N
CDF = CDF + dExpMean * dK
Loop
PoissonInv = Array(N, CDFOld, CDF)
Else
' Plan B, for large means; approximate the Poisson as a normal distribution
' http://en.wikipedia.org/wiki/Continuity_correction#Poisson
Dim iInv As Long
With WorksheetFunction
iInv = .Ceiling(.Norm_Inv(Prob, Mean, Sqr(Mean)) - 0.5, 1)
PoissonInv = Array(iInv, _
.Norm_Dist(iInv - 0.5, Mean, Sqr(Mean), True), _
.Norm_Dist(iInv + 0.5, Mean, Sqr(Mean), True))
End With
End If
End Function
