Calculating miniscule numbers for chi-squared distribution -- numerical precision

Question

I am using the pchisq function in R to calculate the cumulative distribution function for the chi-squared distribution. I would like to calculate very small values, such that 1-pchisq(...) can have a value smaller than 2.2e-16 (which is the numerical precision limit for R's numeric format). Right now, these very small values simply become 0.

Things I've tried:

• Setting the digits display option to 22 (max)

• Using the Rmpfr package for increased precision, but that number format doesn't work with the pchisq function

• Breaking the CDF function into its component gamma functions, but that results in similar precision limits. Any ideas on how I can calculate what I want?

Background: I'm using Fisher's method to combine a bunch of p-values. Yes, I know these p-values are minuscule, but it is actually useful for what I'm analyzing.

Answer 1

A couple of things.

• 2.2e-16 is not the lower precision limit of values in R; it's just the way that R prints very small p-values by default, using format.pval:

format.pval(1e-20)
## [1] "< 2.22e-16"

• values smaller than approximately 1e-320 do round down to zero:

1e-330
## [1] 0

@SeverinPappadeux's suggestion is exactly right:

pchisq(121231,1,lower.tail=FALSE,log.p=TRUE)
## [1] -60621.58

This is equivalent to 10^(-26327):

-60621.58/log(10)
## -26327.62

Check for a less extreme value:

log10(pchisq(100,1,lower.tail=FALSE) )
## [1] -22.81702
pchisq(100,1,lower.tail=FALSE,log.p=TRUE)/log(10)
## [1] -22.81702

Furthermore, log(p) is exactly what you need to use for Fisher's method.