Hi Peter, Rui, Chrstophe and Gabriel, Thanks for your inputs -- the use of qnorm(., log=TRUE) is a good point in line with pnorm with which we devised log(p) as
log(2) + pnorm(-abs(z), lower.tail = TRUE, log.p = TRUE) that could do really really well for large z compared to Rmpfr. Maybe I am asking too much since z <-20000 > Rmpfr::format(2*pnorm(mpfr(-abs(z),100),lower.tail=TRUE,log.p=FALSE)) [1] "1.660579603192917090365313727164e-86858901" already gives a rarely seen small p value. I gather I also need a multiple precision exp() and their sum since exp(z^2/2) is also a Bayes Factor so I get log(x_i )/sum_i log(x_i) instead. To this point, I am obliged to clarify, see https://statgen.github.io/gwas-credible-sets/method/locuszoom-credible-sets.pdf. I agree many feel geneticists go to far with small p values which I would have difficulty to argue againston the other hand it is also expected to see these in a non-genetic context. For instance the Framingham study was established in 1948 just got $34m for six years on phenotypewide association which we would be interesting to see. Best wishes, Jing Hua ________________________________ From: peter dalgaard <pda...@gmail.com> Sent: 21 June 2019 16:24 To: jing hua zhao Cc: Rui Barradas; r-devel@r-project.org Subject: Re: [Rd] Calculation of e^{z^2/2} for a normal deviate z You may want to look into using the log option to qnorm e.g., in round figures: > log(1e-300) [1] -690.7755 > qnorm(-691, log=TRUE) [1] -37.05315 > exp(37^2/2) [1] 1.881797e+297 > exp(-37^2/2) [1] 5.314068e-298 Notice that floating point representation cuts out at 1e+/-308 or so. If you want to go outside that range, you may need explicit manipulation of the log values. qnorm() itself seems quite happy with much smaller values: > qnorm(-5000, log=TRUE) [1] -99.94475 -pd > On 21 Jun 2019, at 17:11 , jing hua zhao <jinghuaz...@hotmail.com> wrote: > > Dear Rui, > > Thanks for your quick reply -- this allows me to see the bottom of this. I > was hoping we could have a handle of those p in genmoics such as 1e-300 or > smaller. > > Best wishes, > > > Jing Hua > > ________________________________ > From: Rui Barradas <ruipbarra...@sapo.pt> > Sent: 21 June 2019 15:03 > To: jing hua zhao; r-devel@r-project.org > Subject: Re: [Rd] Calculation of e^{z^2/2} for a normal deviate z > > Hello, > > Well, try it: > > p <- .Machine$double.eps^seq(0.5, 1, by = 0.05) > z <- qnorm(p/2) > > pnorm(z) > # [1] 7.450581e-09 1.228888e-09 2.026908e-10 3.343152e-11 5.514145e-12 > # [6] 9.094947e-13 1.500107e-13 2.474254e-14 4.080996e-15 6.731134e-16 > #[11] 1.110223e-16 > p/2 > # [1] 7.450581e-09 1.228888e-09 2.026908e-10 3.343152e-11 5.514145e-12 > # [6] 9.094947e-13 1.500107e-13 2.474254e-14 4.080996e-15 6.731134e-16 > #[11] 1.110223e-16 > > exp(z*z/2) > # [1] 9.184907e+06 5.301421e+07 3.073154e+08 1.787931e+09 1.043417e+10 > # [6] 6.105491e+10 3.580873e+11 2.104460e+12 1.239008e+13 7.306423e+13 > #[11] 4.314798e+14 > > > p is the smallest possible such that 1 + p != 1 and I couldn't find > anything to worry about. > > > R version 3.6.0 (2019-04-26) > Platform: x86_64-pc-linux-gnu (64-bit) > Running under: Ubuntu 19.04 > > Matrix products: default > BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.8.0 > LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.8.0 > > locale: > [1] LC_CTYPE=pt_PT.UTF-8 LC_NUMERIC=C > [3] LC_TIME=pt_PT.UTF-8 LC_COLLATE=pt_PT.UTF-8 > [5] LC_MONETARY=pt_PT.UTF-8 LC_MESSAGES=pt_PT.UTF-8 > [7] LC_PAPER=pt_PT.UTF-8 LC_NAME=C > [9] LC_ADDRESS=C LC_TELEPHONE=C > [11] LC_MEASUREMENT=pt_PT.UTF-8 LC_IDENTIFICATION=C > > attached base packages: > [1] stats graphics grDevices utils datasets methods > [7] base > > other attached packages: > > [many packages loaded] > > > Hope this helps, > > Rui Barradas > > �s 15:24 de 21/06/19, jing hua zhao escreveu: >> Dear R-developers, >> >> I am keen to calculate exp(z*z/2) with z=qnorm(p/2) and p is very small. I >> wonder if anyone has experience with this? >> >> Thanks very much in advance, >> >> >> Jing Hua >> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> R-devel@r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-devel >> > > [[alternative HTML version deleted]] > > ______________________________________________ > R-devel@r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel -- Peter Dalgaard, Professor, Center for Statistics, Copenhagen Business School Solbjerg Plads 3, 2000 Frederiksberg, Denmark Phone: (+45)38153501 Office: A 4.23 Email: pd....@cbs.dk Priv: pda...@gmail.com [[alternative HTML version deleted]] ______________________________________________ R-devel@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-devel
