[R] Normal and Poisson tail area expectations in R
Ravi Varadhan
rvaradhan at jhmi.edu
Thu Jun 14 21:04:55 CEST 2007
Inspired by Chuck's elegant solution to the Poisson tail area problem, here
is a simple solution to the normal tail area expectation that does not use
integrate().
Gu.k <- function(k) {1/sqrt(2*pi) * exp(-k*k/2) - k * pnorm(k, lower=FALSE)}
> k <- 1:10
> Gu.k(k)
[1] 8.331547e-02 8.490703e-03 3.821543e-04 7.145258e-06 5.346166e-08
1.563570e-10 1.760326e-13 7.550262e-17 1.224779e-20 7.474560e-25
>
> sapply(k, function(k)integrate(function(x)
(x-k)*dnorm(x),lower=k,upper=Inf)$val)
[1] 8.331547e-02 8.490702e-03 3.821543e-04 7.145259e-06 5.346163e-08
1.563570e-10 1.760326e-13 7.550262e-17 1.224779e-20 7.474560e-25
>
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: rvaradhan at jhmi.edu
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
----------------------------------------------------------------------------
--------
-----Original Message-----
From: r-help-bounces at stat.math.ethz.ch
[mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of Ravi Varadhan
Sent: Thursday, June 14, 2007 2:08 PM
To: 'Charles C. Berry'
Cc: 'kavindra malik'; r-help at stat.math.ethz.ch
Subject: Re: [R] Normal and Poisson tail area expectations in R
Perfect, Chuck. I got a closed-form solution after some algebraic labor,
but your solution is simple and elegant.
Ravi.
----------------------------------------------------------------------------
-------
Ravi Varadhan, Ph.D.
Assistant Professor, The Center on Aging and Health
Division of Geriatric Medicine and Gerontology
Johns Hopkins University
Ph: (410) 502-2619
Fax: (410) 614-9625
Email: rvaradhan at jhmi.edu
Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
----------------------------------------------------------------------------
--------
-----Original Message-----
From: Charles C. Berry [mailto:cberry at tajo.ucsd.edu]
Sent: Thursday, June 14, 2007 1:36 PM
To: Ravi Varadhan
Cc: 'kavindra malik'; r-help at stat.math.ethz.ch
Subject: RE: [R] Normal and Poisson tail area expectations in R
Ravi,
This looks simple to me.
km_G <- function(lambda,k)
lambda*ppois(k-1,lambda,lower=FALSE) -
k*ppois(k,lambda,lower=FALSE)
Am I confused here?
Chuck
On Wed, 13 Jun 2007, Ravi Varadhan wrote:
>
> More interesting is the Poisson convolution. I don't know if there is an
> analytic solution to this. I looked at Jolley's "Summation of Series" and
> Abramowitz and Stegun, but no help there. It seems that discrete FFT
> technique should work. Does anyone know the answer?
>
> Ravi.
>
----------------------------------------------------------------------------
> -------
>
> Ravi Varadhan, Ph.D.
>
> Assistant Professor, The Center on Aging and Health
>
> Division of Geriatric Medicine and Gerontology
>
> Johns Hopkins University
>
> Ph: (410) 502-2619
>
> Fax: (410) 614-9625
>
> Email: rvaradhan at jhmi.edu
>
> Webpage: http://www.jhsph.edu/agingandhealth/People/Faculty/Varadhan.html
>
>
>
>
----------------------------------------------------------------------------
> --------
>
> -----Original Message-----
> From: r-help-bounces at stat.math.ethz.ch
> [mailto:r-help-bounces at stat.math.ethz.ch] On Behalf Of kavindra malik
> Sent: Wednesday, June 13, 2007 5:45 PM
> To: Charles C. Berry
> Cc: r-help at stat.math.ethz.ch
> Subject: Re: [R] Normal and Poisson tail area expectations in R
>
> Thank you very much. This solves the problem I was trying to solve. I am
new
> to R and am learning. A great lesson in the power of R...
>
> "Charles C. Berry" <cberry at tajo.ucsd.edu> wrote: On Wed, 13 Jun 2007,
> kavindra malik wrote:
>
>> I am interested in R functions for the following integrals / sums
> (expressed best I can in text) -
>>
>> Normal: G_u(k) = Integration_{Lower limit=k}^{Upper limit=infinity} [(u
> -k) f(u) d(u)], where where u is N(0,1), and f(u) is the density function.
>>
>> Poisson: G(lambda,k) = Sum_{Lower limit=k}^{Upper limit=infinity} [(x-k)
> p(x, lambda)] where P(x,lambda) is the Poisson prob function with
parameter
> lambda.
>>
>> The Normal expression is very commonly used in inventory management to
>> determine safety stocks (and its tabular values can be found in some
>> texts) - and I am also looking for Poisson and/or Gamma as that'd fit
>> the situation better.
>>
>> I am wondering if there are standard functions in R that might allow me
to
> get these values, instead of needing to do the numerical integration, etc.
> myself.
>
> Not that I know of, but it is not difficult to do the integration:
>
>> k <- 1.1 # for example
>> integrate(function(x) (x-k)*dnorm(x),lower=k,upper=Inf)
> 0.06861951 with absolute error < 5.5e-07
>>
>
> see
>
> ?integrate
> ?qnorm
> ?qpois
> ?qgamma
>
>> Thank you very much.
>>
>>
>>
>> ---------------------------------
>> Sucker-punch spam with award-winning protection.
>>
>> [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at 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
>> and provide commented, minimal, self-contained, reproducible code.
>>
>
> Charles C. Berry (858) 534-2098
> Dept of Family/Preventive
> Medicine
> E mailto:cberry at tajo.ucsd.edu UC San Diego
> http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
>
>
>
>
>
> ---------------------------------
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at 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
> and provide commented, minimal, self-contained, reproducible code.
>
Charles C. Berry (858) 534-2098
Dept of Family/Preventive
Medicine
E mailto:cberry at tajo.ucsd.edu UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901
______________________________________________
R-help at 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
and provide commented, minimal, self-contained, reproducible code.
More information about the R-help
mailing list