[R] Calculate the natural log of cdf between 2 intervals

[Petr Savicky]

On Thu, Feb 02, 2012 at 01:18:42PM -0800, justin jarvis wrote:
> Hello all,
> I was wondering if there is an R function to do the following:
> 
> [*] log(pnorm(x)-pnorm(y)), where x>y.
> 
> I don't want all the area under the natural log of the normal pdf less than
> x, I only want the area between y and x.
> 
> I am aware of the ability to specify log.p=TRUE, which gives me the log of
> the probability that X<=x.  This does not help me, because the following
> code:
> pnorm(x, log.p=TRUE)-pnorm(y,log.p=TRUE) is not the same as [*]
> mathematically.
> 
> I cannot use [*] because some of my x's are far less than the mean, more
> than 10 sd.  This causes me to take the log(0) which is an error.  Thus, I
> need to stay in the log scale, since, for z less than 10 sd below the mean,

Hello:

Try the following.

  x <-  -20
  y <-  -19.9
  xplog <- pnorm(x, log.p=TRUE)
  yplog <- pnorm(y,log.p=TRUE)
  logdiff <- yplog + log1p( - exp(xplog - yplog))
  logdiff

  [1] -202.0626

In an exact arithmetic, we have

  exp(logdiff) =
  exp(yplog + log1p( - exp(xplog - yplog))) =
  exp(yplog + log(1 - exp(xplog - yplog))) =
  exp(yplog) * (1 - exp(xplog - yplog)) =
  exp(yplog) - exp(xplog)

So, we have

  exp(logdiff) = exp(yplog) - exp(xplog)
  logdiff = log(exp(yplog) - exp(xplog))

as required.

Hope this helps.

Petr Savicky.