An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20131130/34ddc066/attachment.pl>
beginner help: 10th index is being overwritten by 11th
3 messages · Loyack, Eric, Berend Hasselman, Duncan Murdoch
On 30-11-2013, at 08:37, Loyack, Eric <eloyack at arcadia.edu> 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
R FAQ 7.31 (http://cran.r-project.org/doc/FAQ/R-FAQ.html#Why-doesn_0027t-R-think-these-numbers-are-equal_003f) Insert print(formatC(10*rho+1,digits=16,format="f")) before the line with print(PTMpvalMeans[rho*10 + 1]) in your code and you should see what?s wrong. Use an explicit counter variable k and increment with k <- k+1 to index your vectors. Berend
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.
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.