[Rd] lower.tail option in pnorm

Ravi Varadhan rvaradhan at jhmi.edu
Tue Dec 8 17:32:13 CET 2009


They will not be the same.  The problem is that the `lower.tail' argument is not vectorized.  Therefore, it is always set equal to the first element of `Resp', which in your example is FALSE.

If you want to obtain same results, this will do the trick:

ans1 <- ifelse(Resp, log(pnorm(Xb)), log(1 - pnorm(Xb)))
ans2 <- apply(cbind(Xb, Resp), 1, function(x) pnorm(x[1], lower.tail=x[2], log.p=TRUE))

all.equal(ans1, as.numeric(ans2))


Ravi.

____________________________________________________________________

Ravi Varadhan, Ph.D.
Assistant Professor,
Division of Geriatric Medicine and Gerontology
School of Medicine
Johns Hopkins University

Ph. (410) 502-2619
email: rvaradhan at jhmi.edu


----- Original Message -----
From: Ken Knoblauch <ken.knoblauch at inserm.fr>
Date: Tuesday, December 8, 2009 10:37 am
Subject: [Rd] lower.tail option in pnorm
To: r-devel at r-project.org


> Hi,
> 
> I would have thought that these two constructions would
> produce the same result but they do not.
> 
> Resp <- rbinom(10, 1, 0.5)
> Stim <- rep(0:1, 5)
> mm <- model.matrix(~ Stim)
> Xb <- mm %*% c(0, 1)
> ifelse(Resp, log(pnorm(Xb)), log(1 - pnorm(Xb)))
> pnorm(as.vector(Xb), lower.tail = Resp, log.p = TRUE)
> 
> > ifelse(Resp, log(pnorm(Xb)), log(1 - pnorm(Xb)))
>   [1] -0.6931472 -1.8410216 -0.6931472 -0.1727538 -0.6931472
>   [6] -0.1727538 -0.6931472 -1.8410216 -0.6931472 -1.8410216
> > pnorm(as.vector(Xb), lower.tail = Resp, log.p = TRUE)
>   [1] -0.6931472 -1.8410216 -0.6931472 -1.8410216 -0.6931472
>   [6] -1.8410216 -0.6931472 -1.8410216 -0.6931472 -1.8410216
> 
> If I have missed something obvious, I would be grateful
> to have it pointed out.
> 
> > sessionInfo()
> R version 2.10.1 beta (2009-12-04 r50668)
> i386-apple-darwin9.8.0
> 
> locale:
> [1] en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8
> 
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods
> [7] base
> 
> loaded via a namespace (and not attached):
> [1] tools_2.10.1
> 
> Thanks for any enlightenment.
> 
> best,
> 
> Ken
> 
> -- 
> Ken Knoblauch
> Inserm U846
> Stem-cell and Brain Research Institute
> Department of Integrative Neurosciences
> 18 avenue du Doyen Lépine
> 69500 Bron
> France
> tel: +33 (0)4 72 91 34 77
> fax: +33 (0)4 72 91 34 61
> portable: +33 (0)6 84 10 64 10
> 
> 
> ______________________________________________
> R-devel at r-project.org mailing list
>



More information about the R-devel mailing list