From: Duncan Murdoch <murdoch_at_stats.uwo.ca>

Date: Thu 21 Jul 2005 - 03:09:52 EST

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Thu Jul 21 03:18:57 2005

Date: Thu 21 Jul 2005 - 03:09:52 EST

On 7/20/2005 12:52 PM, Vincent Goulet wrote:

> We obtained some disturbing results from convolve() (inaccuracies and negative

*> probabilities). We'll try to make the context clear in as few lines as
**> possible...
**>
**> Our function panjer() (code below) basically computes recursively the
**> probability mass function of a compound Poisson distribution. When the
**> Poisson parameter lambda is very large, the starting value of the recursive
**> scheme --- the mass at 0 --- is 0 and the recursion fails. One way to solve
**> this problem is to divide lambda by 2^n, apply the panjer() function and then
**> convolve the result with itself n times.
**>
**> We applied the panjer() function with a lambda such that the mass at 0 is just
**> larger than .Machine$double.xmin. We thus know that once this pmf is
**> convoluted with itself, the first probabilities will be 0 (for the computer).
**>
**> Here are the two issues we have with convolve():
**>
**> 1. The probabilities we know should be 0 are rather in the vicinity of 1E-19,
**> as if convolve() could not "go lower". Using a hand made convolution function
**> (not given here), we obtained the correct values. When probabilities get
**> around 1E-12, results from convolve() and our home made function are
**> essentially identical.
**>
**> 2. We obtained negative probabilities. More accurately, the same example
**> returns negative probabilities under Windows, but not under Linux. We also
**> obtained negative probabilities for another example under Linux, though.
**>
**> We understand that convolve() computes the convolutions using fft(), but we
**> are not familiar enough with the latter to assess if the above issues are
**> some sort of bugs or normal behavior. In the latter case, is there is any
**> workaround?
*

Rounding will depend on the hardware and the compiler, so this might be normal behaviour. It's a little disturbing, but you shouldn't expect calculations to be accurate to more digits than your platform supports. Convolving using an FFT is essentially rotating the vectors in a high dimensional space, multiplying terms, and then rotating back. Unless those rotations are unrealistically accurate very small numbers won't necessarily show up as zeros.

I'd suggest re-thinking your panjer function to work in log probabilities instead of probabilities, so that you can handle a larger dynamic range before you run into underflow problems, and you don't need to use convolve at all. Convert them back to probabilities at the very end.

Duncan Murdoch

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Thu Jul 21 03:18:57 2005

*
This archive was generated by hypermail 2.1.8
: Fri 03 Mar 2006 - 03:33:53 EST
*