[R] beginner help: 10th index is being overwritten by 11th

[Duncan Murdoch]

On 13-11-30 4:37 PM, Loyack, Eric wrote:
> I am using R Studio and writing code to determine p-values for statistical
> t-tests.  Code is below.  When I print out the values in the loop for
> PTMpvalMeans they are correct, but when I store them the 10th one
> overwrites the 11th.  Can someone tell me what I am doing wrong? Thanks in
> advance, Eric

It looks like a fairly obscure case of FAQ 7.31.  It's a bad idea to use 
0.1 as a loop increment, because fractions like that can't be 
represented exactly, and rounding error accumulates.

Duncan Murdoch

>
> rho <- 0.0
> PTMpvalMeans <- numeric(11)
> while (rho < 1.0){
> x1 <- rnorm(25,0,1)
> x2 <- rnorm(25,0,1)
> y1 <- rnorm(25,0,1)
> y2 <- rnorm(25,0,1)
> w1 <- rho*x1 + sqrt(1-rho^2)*y1
> w2 <- rho*x2 + sqrt(1-rho^2)*y2
> controlPretest <- 50 + 10*x1
> groupPretest <- 50 + 10*x2
> controlPosttest <- 60 + 10*w1
> groupPosttest <- 70 + 10*w2
> postTestMethod <-
> t.test(controlPosttest,groupPosttest,mu=0,alternative="two.sided")
> PTMpval <- postTestMethod$p.value
> for (i in 1:999){
> x1 <- rnorm(25,0,1)
> x2 <- rnorm(25,0,1)
> y1 <- rnorm(25,0,1)
> y2 <- rnorm(25,0,1)
> w1 <- rho*x1 + sqrt(1-rho^2)*y1
> w2 <- rho*x2 + sqrt(1-rho^2)*y2
> controlPretestLoop <- 50 + 10*x1
> groupPretestLoop <- 50 + 10*x2
> controlPosttestLoop <- 60 + 10*w1
> groupPosttestLoop <- 70 + 10*w2
> controlPretest <- c(controlPretest,controlPretestLoop)
> groupPosttest <- c(groupPosttest,groupPosttestLoop)
> groupPretest <- c(groupPretest,groupPretestLoop)
> controlPosttest <- c(controlPosttest,controlPosttestLoop)
> postTestMethodLoop <-
> t.test(controlPosttestLoop,groupPosttestLoop,mu=0,alternative="two.sided")
> PTMpvalLoop <- postTestMethodLoop$p.value
> PTMpval <- c(PTMpval,PTMpvalLoop)
> }
> PTMpvalMeans[rho*10 + 1] <- mean(PTMpval)
> print(rho)
> print(PTMpvalMeans[rho*10 + 1])
> rho <- rho + 0.1
> }
> PTMpvalMeans
>
> 	[[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>