Tuesday, February 17, 2015

R22. One Sample t-test in R

A one-sample t-test is conducted over a range of mean values. The null hypothesis is that test mean is the true mean of the population.


The sample is of length 500 and is a random sample based on a standard normal distribution, which has a mean of 0. Thus, we expect values far from 0 will be rejected.


If you uncomment the cat() lines, you will get intermediate results. You can use the t.test function or a calculation, based on the t-distribution. We have a two sided alternative hypotheses, that is, we sum over both tails. This is the default for a t.test.

# ex22.R
N <- 500
mu <- seq(-.2,+.2,.001) # test for the population mean
pr <- numeric(length(mu))
rej <- numeric(length(mu))
samp <- rnorm(N)
mean.samp <- mean(samp)
sd.samp <- sd(samp)
for (i in 1:length(mu)) {
  m <- mu[i]
  t.crit <- (mean.samp-m)/(sd.samp/sqrt(N))
  if (t.crit>=0) {
    p <- 2*(1-pt(t.crit, df = N-1))
  } else {
    p <- 2*pt(t.crit, df = N-1 )
  }
  pr[i] <- p
  res <- t.test(samp, mu = m)
  #cat('\n-- t.crit = ',t.crit,'\tres$statistic =',res$statistic)
  #cat('\n-- p =',p,'\tres$p.value =', res$p.value)
  if (p<0.05) {
    #cat('\nmu =',m,' rejected')
    rej[i] <- 1
  } else {
    #cat('\nmu =',m,' not rejected')
    rej[i] <- 0
  }
}
plot(mu, rej, type = 'l',
     main = "Means rejected (p-value: blue), (threshold: green)",
     col = 'red', xlab = 'mean', ylab = "Reject = 1 (red)")
lines(mu, pr, col = 'blue')
abline(h = .05, col = 'green')

Output:

No comments:

Post a Comment