(35S) Statistical Distributions Functions
11-21-2015, 07:27 PM (This post was last modified: 11-27-2015 06:46 AM by Dieter.)
Post: #74
 Dieter Senior Member Posts: 2,398 Joined: Dec 2013
RE: HP 35s Statistical Distributions Functions
(11-20-2015 11:06 PM)Dieter Wrote:  I have reworked this a bit so that now the same values as with the other distributions are returned. This includes the upper CDF (which in these cases is not just 1 minus the lower CDF) as well as the PDF (or PMF - probability mass function -, to be more precise).
(...)
Edit: it looks like the Poisson quantile is not that complicated. At least an experimental version seems to run without major problems. ;-)

OK, here is a first version of the complete Poisson distribution: upper and lower CDF, PMF and – the inverse cdf (quantile). Usage is the same as with other programs of this pack. Only the quantile function returns something different:

As the Poisson distribution handles a discrete random variable (x=0, 1, 2, 3, ...) the returned quantile is an integer. However, in most cases there is no exact integer for a given cumulative probability. That's why usually the quantile is defined as the smallest integer with a cdf greater than or equal to the given probability. For instance, if p=0,9 while cdf(11)=0,85 and cdf(12)=0,93 the quantile is 12.

This may lead to off-by-one errors in very close cases where the calculated cdf may be not exact in the last digits. That's why the following program returns three values. In x the calculated quantile is returned, and y at the same time holds the corresponding probability (cdf). Scrolling down the stack with R↓ displays the cdf for the next lower integer, i.e. for x–1. This way you can see the cdfs for the two adjacent integers that bracket the given probability.

Example: mean = 7,3 and p=0,9. The program returns x=11 as well as 0,9319 = cdf(11) and 0,8788 = cdf(10).

The method is quite simple and straightforward. It implements some insights of a 2013 paper by Michael Short titled "Improved Inequalities for the Poisson and Binomial Distribution and Upper Tail Quantile Functions". Actually the Poisson quantile can (almost) be calculated directly by means of Lambert's W function, but since the 35s does not offer one (unlike the 34s) I chose a different approach.

Finally, here is the code:

Code:
P001  LBL P     Poisson distribution starts here P002  XEQ J012  initialize P003  INPUT M   prompt for mean P004  X<=0? P005  GTO P003  reject non-positive entry P006  XEQ A001  prompt for x resp. p P007  FS? 1     quantile mode? P008  GTO P040  go to quantile function P009  X<0? P010  GTO P006  reject negative x P011  RCL+ O P012  STO A     a = x+1 P013  RCL M P014  STO Y     y = mean P015  XEQ G006  call CDF = incomplete Gamma(x+1, mean) P016  XEQ P031  call PMF(x) P017  STO C P018  +         add PMF to (1-CDF) to get upper CDF P019  x<> C     store in C P020  STO D     and move PMF to D P021  x<>y P022  RCL C     arrange data for output P023  x<>y P024  ENTER P025  CLx P026  STO A     clear A and B P027  STO B P028  R↓        and T as well P029  STOP      return T: 0  Z: PMF  Y: upper CDF  X: lower CDF P030  GTO P006 P031  RCL M     calculate PMF P032  LN        use logs to avoid overflow P033  RCLx X P034  RCL- M P035  e^x P036  RCL X P037  !         ln Gamma would be nice here P038  ÷         = PMF P039  RTN P040  RCL M     start of quantile function P041  +/- P042  e^x P043  RCL P P044  x>y?      p > CDF(0) ? P045  GTO P050  then continue P046  CLx P047  STO X     else return 0 P048  x<>y P049  GTO P104  jump to exit P050  XEQ E001  get Normal quantile estimate P051  RCLx E P052  3 E-4     make sure Normal quantile P053  +         is slightly high P054  ENTER P055  x² P056  6 P057  ÷ P058  x<>y P059  RCL M P060  √x P061  x P062  + P063  RCL+ M P064  IP P065  x<=0? P066  RCL T P067  RCL+ O P068  STO X     first guess, may be slightly high P069  RCL+ O P070  STO A P071  RCL M P072  STO Y P073  XEQ G006  calculate CDF P074  x>y?      CDF < 0,5? P075  GTO P080 P076  +/-       if not, p is closer to 1 P077  RCL P     so better compare -(1-CDF) P078  RCL- O    and -(1-p) P079  GTO P082 P080  R↓ P081  RCL P P082  x<>y P083  STO W P084  XEQ P031  calculate PMF(x) P085  -         CDF(x) - PMF(x) = CDF(x-1) P086  x<y?      is this less than p ? P087  GTO P097  then exit P088  STO W     else save CDF(x) P089  LASTx P090  RCLx X    calculate PMF(x) from PMF(x+1) P091  RCL÷ M P092  DSE X     decrement x P093  GTO P085  and continue P094  -         if ever executed, this should round to zero P095  CLx       so make it exactly zero ;-) P096  ENTER P097  RCL W     now x and y hold the previous and current CDF P098  x>0?      are we working with the CDFs (not with 1-CDF) ? P099  GTO P104 P100  RCL+ O    if not, restore CDF = 1+[-(1-CDF)] P101  x<>y P102  RCL+ O    do so both for CDF(x) and CDF(x-1) P103  x<>y P104  RCL X P105  ENTER P106  CLx P107  R↓        clear T P108  STOP      return T: 0  Z: CDF(x-1)  Y:CDF(x)  X: quantile x P109  GTO P006

Try it and see what you get. As usual, this is supposed to work with the functions J, A, E and G included in the distributions package that was posted as an .ep file for the 35s emulator. I did not do much testing, so all error reports are welcome.

Edit (2015-11-23): added line P065/P066. However, there are still a few cases with larger mean an very small p that may return a slightly low result, e.g. m=50 and p=1E-17 returns 2 instead of 4.

Edit2 (2015-11-26): this problem can be avoided by a different approach which also seems to work for the Binomial quantile.

Dieter
 « Next Oldest | Next Newest »

 Messages In This Thread (35S) Statistical Distributions Functions - PedroLeiva - 10-25-2015, 12:19 AM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-25-2015, 01:29 PM RE: HP 35s Statistical Distributions Functions - SlideRule - 10-26-2015, 03:03 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-26-2015, 04:41 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-26-2015, 08:30 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-26-2015, 10:29 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-26-2015, 11:04 PM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-26-2015, 08:39 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-26-2015, 09:25 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-27-2015, 02:29 PM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-26-2015, 09:28 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-27-2015, 11:59 AM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-27-2015, 04:48 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-27-2015, 05:22 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-27-2015, 07:32 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-27-2015, 10:25 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 12:24 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-28-2015, 12:34 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 06:59 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-27-2015, 11:48 PM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-27-2015, 10:29 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-27-2015, 11:41 PM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-28-2015, 12:11 AM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 12:36 AM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-28-2015, 01:05 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-28-2015, 11:55 AM RE: HP 35s Statistical Distributions Functions - Thomas Klemm - 10-28-2015, 01:40 AM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 06:21 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-28-2015, 07:54 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 10:34 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-30-2015, 02:36 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-30-2015, 01:56 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-30-2015, 05:20 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-30-2015, 07:16 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-30-2015, 07:43 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-30-2015, 11:08 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-30-2015, 11:16 PM RE: HP 35s Statistical Distributions Functions - rprosperi - 10-31-2015, 09:11 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-28-2015, 10:45 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 10-29-2015, 12:09 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-29-2015, 07:01 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-29-2015, 07:18 AM RE: HP 35s Statistical Distributions Functions - Dieter - 10-29-2015, 11:56 PM RE: HP 35s Statistical Distributions Functions - Dieter - 10-31-2015, 09:05 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-01-2015, 12:44 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-01-2015, 09:17 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-05-2015, 01:37 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-05-2015, 02:16 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-13-2015, 06:59 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-13-2015, 08:04 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-13-2015, 09:05 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-14-2015, 11:22 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-15-2015, 03:04 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-15-2015, 09:26 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-16-2015, 01:13 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-16-2015, 07:08 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-16-2015, 10:47 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-17-2015, 12:44 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-17-2015, 09:20 PM RE: HP 35s Statistical Distributions Functions - Gerald H - 11-19-2015, 09:24 AM RE: HP 35s Statistical Distributions Functions - Dieter - 11-19-2015, 01:00 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-19-2015, 01:17 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-19-2015, 08:31 PM RE: HP 35s Statistical Distributions Functions - Gerald H - 11-19-2015, 05:52 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-20-2015, 12:09 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-20-2015, 01:19 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-20-2015, 05:45 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-20-2015, 08:30 PM RE: HP 35s Statistical Distributions Functions - Marcio - 11-20-2015, 08:37 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-20-2015, 08:54 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-20-2015, 11:06 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-21-2015 07:27 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-20-2015, 05:47 PM RE: HP 35s Statistical Distributions Functions - Marcio - 11-20-2015, 12:24 PM RE: HP 35s Statistical Distributions Functions - PedroLeiva - 11-26-2015, 11:07 PM RE: HP 35s Statistical Distributions Functions - Dieter - 11-27-2015, 06:41 AM RE: HP 35s Statistical Distributions Functions - Dieter - 11-29-2015, 06:50 PM

User(s) browsing this thread: 1 Guest(s)