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 <jinghuazhao at 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 <ruipbarradas at sapo.pt> > Sent: 21 June 2019 15:03 > To: jing hua zhao; r-devel at 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 at r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-devel >> > > [[alternative HTML version deleted]] > > ______________________________________________ > R-devel at 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.mes at cbs.dk Priv: PDalgd at gmail.com
Hi Jing, Peter pointed out how you can, more or less, get numbers for this, and he's absolutely right. At the risk of giving unsolicited advice, though, Im don't think you *should* in this case. Someone on this list with more applied Statistics or Statistical Genetics experience can correct me if I'm wrong, but I have always been under the impression that "pvalues" that small (and honestly ones that are not nearly that small) are not meaningful as actual numbers and should never be used as such. Essentially, AFAIK, a p-value of 1e-300 just means "really small"/"really unlikely", the same as a pvalue of 1e-280 does (even though the second pvalue is *20 orders of magnitude* larger). Like, how *precicely* correct would the distributional assumptions the test is making have to be for that -300 to be meaningful? Note that if you have putative pvalues in the 1e-300 range, a full order of magnitude increase in the pvalue (to 1e-299) only corresponds to only a reduction of 0.06 in the quantile value. Is that amount of measurement difference likely to be meaningful in a genomics setting or is it probably machine noise? (my impression is the later).> log(1e-300)[1] -690.7755> log(1e-299)[1] -688.4729> qnorm(-691, log = TRUE)[1] -37.05315> qnorm(-688, log=TRUE)[1] -36.97216> qnorm(-688.4729, log=TRUE) - qnorm(-690.7755, log = TRUE)[1] 0.06216021 And all of the above ignores the most important thing in practice, though (sorry for burying the lede), which is the fact that the precision tolerance of the computation calculating the pvalue is very likely *hundreds of orders of magnitude *larger (ie less precise) than 1e-300 unless it uses infinite precision arithmetic (which, I think, would be a very weird thing to do for the reasons pointed out above). Personally, whenever I See values on the order of 1e-K where K is bigger than like 20 on the ouside , I just see rounding errors around 0, not numbers that are, e.g., meaningfully comparable or sortable or summable, etc. So R will give you numbers for this, as Peter showed you, but personally remain pretty skeptical that they will actually be useful for what you want to do with them i this case. Best, ~G On Fri, Jun 21, 2019 at 9:32 AM peter dalgaard <pdalgd at gmail.com> wrote:> 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 <jinghuazhao at 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 <ruipbarradas at sapo.pt> > > Sent: 21 June 2019 15:03 > > To: jing hua zhao; r-devel at 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 at r-project.org mailing list > >> https://stat.ethz.ch/mailman/listinfo/r-devel > >> > > > > [[alternative HTML version deleted]] > > > > ______________________________________________ > > R-devel at 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.mes at cbs.dk Priv: PDalgd at gmail.com > > ______________________________________________ > R-devel at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel >[[alternative HTML version deleted]]
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 <pdalgd at gmail.com> Sent: 21 June 2019 16:24 To: jing hua zhao Cc: Rui Barradas; r-devel at 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 <jinghuazhao at 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 <ruipbarradas at sapo.pt> > Sent: 21 June 2019 15:03 > To: jing hua zhao; r-devel at 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 at r-project.org mailing list >> https://stat.ethz.ch/mailman/listinfo/r-devel >> > > [[alternative HTML version deleted]] > > ______________________________________________ > R-devel at 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.mes at cbs.dk Priv: PDalgd at gmail.com [[alternative HTML version deleted]]
I agree with many the sentiments about the wisdom of computing very
small p-values (although the example below may win some kind of a prize:
I've seen people talking about p-values of the order of 10^(-2000), but
never 10^(-(10^8)) !). That said, there are a several tricks for
getting more reasonable sums of very small probabilities. The first is
to scale the p-values by dividing the *largest* of the probabilities,
then do the (p/sum(p)) computation, then multiply the result (I'm sure
this is described/documented somewhere). More generally, there are
methods for computing sums on the log scale, e.g.
https://docs.scipy.org/doc/scipy-0.14.0/reference/generated/scipy.misc.logsumexp.html
I don't know where this has been implemented in the R ecosystem, but
this sort of computation is the basis of the "Brobdingnag" package for
operating on very large ("Brobdingnagian") and very small
("Lilliputian") numbers.
On 2019-06-21 6:58 p.m., jing hua zhao wrote:> 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 <pdalgd at gmail.com>
> Sent: 21 June 2019 16:24
> To: jing hua zhao
> Cc: Rui Barradas; r-devel at 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 <jinghuazhao at
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 <ruipbarradas at sapo.pt>
>> Sent: 21 June 2019 15:03
>> To: jing hua zhao; r-devel at 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 at r-project.org mailing list
>>> https://stat.ethz.ch/mailman/listinfo/r-devel
>>>
>>
>> [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-devel at 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.mes at cbs.dk Priv: PDalgd at gmail.com
>
>
>
>
>
>
>
>
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-devel at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-devel
>
On 22/06/2019 00:58, jing hua zhao wrote:> Hi Peter, Rui, Chrstophe and Gabriel, > > Thanks for your inputs -- the use of qnorm(., log=TRUE) is a good pointAnother approach could be simply to note that a function defined as f(p)=exp(-z(p)^2/2) is regular around p=0 with f(0)=0. It has roughly the shape of p*(2-p) for p \in [0; 1]. So we can calculate let say f(10^-10) with sufficient precision using Rmpfr and then use a linear approximation for p from [0, 10^-10]. After that a simple inverse gives us e^(z*z/2). Serguei.> 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 <pdalgd at gmail.com> > Sent: 21 June 2019 16:24 > To: jing hua zhao > Cc: Rui Barradas; r-devel at 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 <jinghuazhao at 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 <ruipbarradas at sapo.pt> >> Sent: 21 June 2019 15:03 >> To: jing hua zhao; r-devel at 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 at r-project.org mailing list >>> https://stat.ethz.ch/mailman/listinfo/r-devel >>> >> [[alternative HTML version deleted]] >> >> ______________________________________________ >> R-devel at 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.mes at cbs.dk Priv: PDalgd at gmail.com > > > > > > > > > > > [[alternative HTML version deleted]] > > ______________________________________________ > R-devel at r-project.org mailing list > https://stat.ethz.ch/mailman/listinfo/r-devel