In [1]:

    

boottest1 <- function(x, theta0, B) {
  n <- length(x)
  v <- mean(x)
  z <- x - v + theta0
  counter <- 0
  teststatall <- rep(0,B)
  for (i in 1:B) {
    xstar <- sample(z, n, replace=TRUE)
    vstar <- mean(xstar)
    if (vstar >= v) counter <- counter + 1
    teststatall[i] <- vstar
  }
  pvalue <- counter/B
  #list(origtest=v, pvalue=pvalue, teststatall=teststatall)
  list(origtest=v, pvalue=pvalue)
}


    

In [2]:

    

boottest1.twosided <- function(x, theta0, B) {
  n <- length(x)
  v <- mean(x)
  z <- x - v + theta0
  counter <- 0
  teststatall <- rep(0,B)
  for (i in 1:B) {
    xstar <- sample(z, n, replace=TRUE)
    vstar <- mean(xstar)
    if (abs(vstar-theta0) >= abs(v-theta0)) counter <- counter + 1
    teststatall[i] <- vstar
  }
  pvalue <- counter/B
  #list(origtest=v, pvalue=pvalue, teststatall=teststatall)
  list(origtest=v, pvalue=pvalue)
}


    

In [3]:

    

rcn <- function(n, eps, sigmac) {
  ind <- rbinom(n, 1, eps)
  x <- rnorm(n)
  rcn <- x*(1-ind)+sigmac*x*ind
  rcn
}


    

In [4]:

    

w <- rcn(20, 0.2, 4)
x <- 10 * w + 100


    

In [5]:

    

boottest1(x, 90, 3000)


    

    
Out[5]:





	
$origtest
		
102.010082502582
	
$pvalue
		
0.0546666666666667

In [6]:

    

t.test(x, mu=90, alternative="greater")


    

    
Out[6]:





	One Sample t-test

data:  x
t = 1.6614, df = 19, p-value = 0.05652
alternative hypothesis: true mean is greater than 90
95 percent confidence interval:
 89.51069      Inf
sample estimates:
mean of x 
 102.0101 

In [7]:

    

boottest1.twosided(x, 90, 3000)


    

    
Out[7]:





	
$origtest
		
102.010082502582
	
$pvalue
		
0.0786666666666667

In [8]:

    

boottest1.med <- function(x, theta0, B) {
  n <- length(x)
  v <- median(x)
  z <- x - v + theta0
  counter <- 0
  teststatall <- rep(0,B)
  for (i in 1:B) {
    xstar <- sample(z, n, replace=TRUE)
    vstar <- median(xstar)
    if (vstar >= v) counter <- counter + 1
    teststatall[i] <- vstar
  }
  pvalue <- counter/B
  #list(origtest=v, pvalue=pvalue, teststatall=teststatall)
  list(origtest=v, pvalue=pvalue)
}

boottest1.med.twosided <- function(x, theta0, B) {
  n <- length(x)
  v <- median(x)
  z <- x - v + theta0
  counter <- 0
  teststatall <- rep(0,B)
  for (i in 1:B) {
    xstar <- sample(z, n, replace=TRUE)
    vstar <- median(xstar)
    if (abs(vstar-theta0) >= abs(v-theta0)) counter <- counter + 1
    teststatall[i] <- vstar
  }
  pvalue <- counter/B
  #list(origtest=v, pvalue=pvalue, teststatall=teststatall)
  list(origtest=v, pvalue=pvalue)
}


    

In [9]:

    

boottest1.med(x, 90, 3000)


    

    
Out[9]:





	
$origtest
		
95.8467802474429
	
$pvalue
		
0.146333333333333

In [10]:

    

boottest1.med.twosided(x, 90, 3000)


    

    
Out[10]:





	
$origtest
		
95.8467802474429
	
$pvalue
		
0.161333333333333

In [ ]: