# 
# Written by:
# -- 
# John L. Weatherwax                2005-08-04
# 
# email: wax@alum.mit.edu
# 
# Please send comments and especially bug reports to the
# above email address.
# 
# Epage 169
# 
#-----

X1 = c(3.5, 3.7, 3.9, 4.1, 3.4, 3.5, 4.1, 3.8, 3.6, 3.7)
X2 = c(3.2, 3.6, 3.1, 3.4, 3.0, 3.4, 2.8, 3.1, 3.3, 3.6)

xbar = mean(X1)
sigma1Hat2 = var(X1) 
ybar = mean(X2)
sigma2Hat2 = var(X2)

N = length(X1)
sz2 = 0.5 * ( sigma1Hat2 + sigma2Hat2 )

q = ( xbar-ybar ) / ( sqrt(sz2)*sqrt(2/N) )

rho = 0.05
ll = qt( 0.5*rho, 2*N-2 )
ul = qt( 1-0.5*rho, 2*N-2 )

print( sprintf("q= %10.6f", q) )
print( sprintf("acceptance interval (for H_0)= [%10.6f,%10.6f]",ll,ul) )