[WP34s] Inverse discrete probability distributions  Printable Version + HP Forums (https://archived.hpcalc.org/museumforum) + Forum: HP Museum Forums (https://archived.hpcalc.org/museumforum/forum1.html) + Forum: Old HP Forum Archives (https://archived.hpcalc.org/museumforum/forum2.html) + Thread: [WP34s] Inverse discrete probability distributions (/thread221599.html) 
[WP34s] Inverse discrete probability distributions  Eduardo Duenez  05152012 Hello all, The post will be long, but I see no other way to make what I think is a very important point clearly enough. I was helping my daughter with her homework. She had to spin a roulette many times and count the frequencies of different outcomes. I wrote a 34s program to do a Monte Carlo simulation. That way she could do about 1000 rolls in 2 seconds (a loop calling RAN# and increasing the appropriate register each time). Just for fun, I decided to write a second program to do the same but using the inverse binomial distribution instead, only to quickly realize that the implementation of inverse discrete distributions on the 34s is "off by one". Let me explain. For simplicity let us consider the geometric distribution, say with parameter p=0.5 stored in J. Such a geometric variable takes the value m=0 with probability 1/2, m=1 with probability 1/4, m=2 with prob=1/8, etc. The distribution function (the function Geom in the 34s), is the function F(x) defined for all real values of x as follows:
F(x) = 0 if x<0 The "inverse distribution" F^(1)(x), that is, the inverse function to F(x) above is, strictly speaking only defined for all real values of x, but only for the values of x=0,1/2,3/4,7/8,..., etc. And worse, it is not *uniquely* defined for those values. For instance, we could take F^(1)(1/2) to be equal to any number x such that 0<=x<1. This is mathematically correct, strictly speaking. The question becomes: What is the more *useful* way to define the inverse functions? The 34s calls that inverse distribution function above Geom^{1}. Currently, the 34s calculates this function as follows:
Geom^{1}(x) = 0 if x<0, !!! The output for x<0 notwithstanding (which I propose to make simply undefined), the values returned for x>=0 are almost always one unit smaller than they ought to be. The values returned should be:
Geom^{1}(x) is undefined if x<0 (one can probably argue this should be either 1, or infinity, or NaN. The point is that nobody should ever need to evaluate an inverse distribution function when the argument is not in the range 0<=x<=1.), Unfortunately, it seems as if all of the discrete inverse distributions are "off by 1" in the 43s currently (Caveat: I only checked the geometric and binomial). I want to point out an important immediate corollary of fixing the inverse distributions as I propose above: One can use the inverse distribution (that is, the "inverse transform method" as explained in these notes) to go from a uniform RAN# variable to a variable with the desired distribution:
RAN# generates a random variable distributed geometrically. This method fails with the current implementation of the inverse discrete random variables. If we replace Geom^{1} by, say, Binom^{1}, Expon^{1} or Phi^{1} we get a binomially, exponentially or normally distributed random variable, regardless of whether the distribution is discrete or continuous. Elegant and consistent. Eduardo PS: I have the source code of the probability distributions, but I could not quickly make sense of it so as to point out more specifically how to implement the changes. Sorry. [Edited to exchange some <= to < in order that the examples I give indeed be inverse to the original distribution.]
Edited: 15 May 2012, 11:40 p.m. after one or more responses were posted
Re: [WP34s] Inverse discrete probability distributions  Walter B  05152012 Eduardo,
Thanks for the link and for reporting your observations. It will take us some time to check and test, but we will return to you as soon as we know more :)
Re: [WP34s] Inverse discrete probability distributions  Eduardo Duenez  05152012 Actually, the link I provided is probably not a very good reference on the "inverse distribution method", but merely the first document a quick search turned up. However, this is a standard result found in just about any textbook on prob/stats. If you want me to send a more detailed/technical explanation (or a better reference) I can provide it, though perhaps the forum may not be the best place for the discussion of the details. (Can always reach me by email: firstname dot lastname at gmail dot com)
Eduardo
Re: [WP34s] Inverse discrete probability distributions  Paul Dale  05152012 This is by design. The inverse of the discrete functions returns, for a given probability p and distribution X, the largest integer n such that the Pr(X <= n) <= p. I.e. it is rounding the value associated with the probability down. Thus for 0 which has Pr(X=0) of 0.5, the inverse won't return 0 until p >= 0.5. What you seem to want is for the inverse to round up. I.e. the smallest integer n such that Pr(x <= n) > p. Both approaches are valid, the question is which one is more desirable here. Anyone want to make a case either way? The x<0 case is a bug and it will be fixed once this question is answered.  Pauli
Edited: 15 May 2012, 9:03 p.m.
Re: [WP34s] Inverse discrete probability distributions  Paul Dale  05152012 Quote: The final fix for the discrete distributions is at label qf_discrete_final in trunk/xrom/distributions/geometric.wp34s. Changing this is the easiest way to alter the rounding behaviour of all discrete distributions at once.
Re: [WP34s] Inverse discrete probability distributions  Eduardo Duenez  05152012 Well, my original post advocates the "rounding up" method, and my opinion is strongly for it. Although I'm no Wolfram fanboy, they have a decent discussion on the "method of inverse transforms" http://demonstrations.wolfram.com/TheMethodOfInverseTransforms/ The "Details" section of that page explains how to use the inverse distribution function to generate a variable with any given distribution starting from a uniformly distributed variable as simulated by RAN#. I want to emphasize that this method *only* works if the inverse distribution is computed rounding up rather than down. I consider the use of the inverse distribution function in transforming uniformly distributed RAN# numbers into random numbers with a *different* distribution to be the quintessential use of these inverses (when "rounding up" is used, as I have advocated, in the way I find mathematically more natural). So I ask the question: Under what circumstances would it be more natural or useful to define the inverse distribution by rounding down? If I have a chance tomorrow, I'll try to find a proper reference (e.g., a textbook) on this matter.
Eduardo
Re: [WP34s] Inverse discrete probability distributions  Eduardo Duenez  05162012 For reference, here's a couple of web pages describing how the inverse distribution function should be defined (either of these is far better than the reference I gave in my original post, which although correct is just lecture notes for a sophomorelevel course in actuarial models). First one, by faculty at the Computer Science Dept. of Ben Gurion University in Israel: http://www.cs.bgu.ac.il/~mps042/invtransnote.htm Another one, written in complete generality (using slightly higher mathematical lingo and including a proof of the Inverse Distribution Theorem), by Ursa Pantle, a statistician at Universität Ulm in Germany: http://www.mathematik.uniulm.de/stochastik/lehre/ss06/markov/skript_engl/node29.html In a nutshell, for a variable with distribution F(x) = probability that the random variable takes values <=x, the generalized inverse distribution function is defined by: F^{1}(y) = infimum of values x such that F(x)>=y. In particular, if F(x) < y then F^{1}(y) must be at least x, and in fact if the underlying distribution is discrete, then F^{1}(y)>x in other words, rounding must be done up, never down.
Eduardo
Re: [WP34s] Inverse discrete probability distributions  Paul Dale  05162012 These links are talking about generating distributed random numbers  this doesn't necessary relate to the inverse cumulative distribution function for a discrete random variable. The inverse function is ill defined  both what the 34S and what you want are inverses in the sense that they return the correct values over the integers which is the domain of the discrete distributions. That all said, nobody else has said anything thus far and I really don't mind which way these functions are defined (it is only an INC or DEC between them), so I'll go ahead and change them to the upward rounding flavours in due course.
Re: [WP34s] Inverse discrete probability distributions  Eduardo Duenez  05162012 Quote: Precisely, since the strict inverses are undefined for most real numbers, it becomes a matter of extending the definition in the most useful manner. Probabilists and statisticians have a preferred way of defining this generalized inverse (viz., the documents I have referenced in my prior posts in this thread). I noticed the 34s inverse discrete distribution functions differ from such definition in most of their domain, and I have tried pointing out so. That said, I certainly do not wish to impose the work on you or the 34s team. I'm sure I can figure out how to make the code changes myself and contribute them to the project given a little time. So long as the changes are welcome, that is.
Eduardo
Re: [WP34s] Inverse discrete probability distributions  Paul Dale  05162012 The changes are already done, just not built yet. They'll need testing however... The only constraint I remember on discrete distributions was they had to be right continuous and both of these definitions are that. Well, beyond always being nonnegative and integrating/summing to unity as all distributions must. It's been a while since I classed myself a statistician and longer still since I studied postgraduate probability theory so I might have forgotten something....
Re: [WP34s] Inverse discrete probability distributions  Eduardo Duenez  05172012 Quote: Consider me first in line to test! (In due course.)
Eduardo
Re: [WP34s] Inverse discrete probability distributions  fhub  05172012 Quote:Well, you can test! ;)
The changes are already builtin:
Franz
