From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: (qmail 1760 invoked by alias); 6 Jul 2007 09:31:01 -0000 Received: (qmail 1747 invoked by uid 22791); 6 Jul 2007 09:30:57 -0000 X-Spam-Check-By: sourceware.org Received: from mail.cs.tut.fi (HELO mail.cs.tut.fi) (130.230.4.42) by sourceware.org (qpsmtpd/0.31) with ESMTP; Fri, 06 Jul 2007 09:30:51 +0000 Received: from spam2.cs.tut.fi (spam2.cs.tut.fi [130.230.4.7]) by mail.cs.tut.fi (Postfix) with ESMTP id 6586A5525; Fri, 6 Jul 2007 12:30:06 +0300 (EEST) Received: from mail.cs.tut.fi ([130.230.4.42]) by spam2.cs.tut.fi (spam2.cs.tut.fi [130.230.4.7]) (amavisd-maia, port 10024) with ESMTP id 13393-22-7; Fri, 6 Jul 2007 12:30:05 +0300 (EEST) Received: from modeemi.modeemi.cs.tut.fi (modeemi.modeemi.cs.tut.fi [130.230.72.134]) by mail.cs.tut.fi (Postfix) with ESMTP id BA7E91282; Fri, 6 Jul 2007 12:30:05 +0300 (EEST) Received: from jolt.modeemi.cs.tut.fi (jolt.modeemi.cs.tut.fi [130.230.72.144]) by modeemi.modeemi.cs.tut.fi (Postfix) with ESMTP id 9DD222C994; Fri, 6 Jul 2007 12:30:05 +0300 (EEST) Received: by jolt.modeemi.cs.tut.fi (Postfix, from userid 16311) id D2A4D4F940; Fri, 6 Jul 2007 12:30:03 +0300 (EEST) Date: Fri, 06 Jul 2007 09:31:00 -0000 To: Jason Detwiler Cc: gsl-discuss@sources.redhat.com Subject: Re: request for Poisson CDF inverse Message-ID: <20070706093002.GI24297@jolt.modeemi.cs.tut.fi> References: MIME-Version: 1.0 Content-Type: text/plain; charset=iso-8859-1 Content-Disposition: inline In-Reply-To: User-Agent: Mutt/1.5.11 From: shd@jolt.modeemi.cs.tut.fi (Heikki Orsila) Mailing-List: contact gsl-discuss-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Subscribe: List-Archive: List-Post: List-Help: , Sender: gsl-discuss-owner@sourceware.org X-SW-Source: 2007-q3/txt/msg00005.txt.bz2 On Thu, Jul 05, 2007 at 05:16:38PM -0700, Jason Detwiler wrote: > I am looking to calculate the inverse of the CDF for the Poisson > distribution. I noticed that many of the CDFs coded in gsl have > inverse functions also coded, but this is not the case for the Poisson > distribution (unfortunately for me). I would like to "request" these > functions. I would try to contribute them myself, but I'm not > confident I could code them to the standards of gsl. What are the > chances of my wishes being fulfilled by another developer? Are you looking for something like this? http://zakalwe.fi/~shd/code/poisson.py This code maps a probability from range [0, 1) into a number of random events happening at a given invocation. This is done in Poisson_Process.random(). Example, 1E6 invocations, probability p = 0.05 of an event. This will give frequencies for each N (the number of random events happening): $ ./poisson.py |sort -n |uniq -c |sort -n ... 75 3 2482 2 47754 1 949688 0 As a comparison: 1E6*0.05^1 = 50000 1E6*0.05^2 = 2500 1E6*0.05^3 = 125 So these values are pretty close to what they should be. This random process is only an approximation. I put an intentional limit into probabilities so that events that rarer than 1E-9 do not happen unless the Poisson_Process constructor is suplied with a value. -- Heikki Orsila Barbie's law: heikki.orsila@iki.fi "Math is hard, let's go shopping!" http://www.iki.fi/shd