Question

Question 6 Consider two independent samples from 3-variate multivariate normal populations: 3 2 2 S2 3 2 1411 Population 1 with 12: sample size n 10, 13 3 2 21 Population 2 with 2 22sample size n2 - 10, 6,S2 2 3 2 μ21 23 Furthermore, assume the population covariance matrices of the two populations are the same. We aim to test 1. Test Ho with a 0.05 by Hotelling's T2 2. Test Ho with α 0.05 by Bonferroni correction. Hint: Find C such that the null can be written as Ho Ci C2

Answer 1

eaY-N_n (ly.nt), c , C2c'たと.79 ว่า

matthod is qvenas 2.0 1.8 o ② C-1.80,-0.20)

Corresponding R code is given as;

x1_bar=c(3,4,5)
x2_bar=c(4,6,6)

S1=matrix(c(3,2,2,2,3,2,2,2,2),ncol=3,byrow=TRUE)
S2=matrix(c(3,2,2,2,3,2,2,2,4),ncol=3,byrow=TRUE)

C=matrix(c(1,-1,0,0,1,-1,1,0,-1),nrow=3,byrow=TRUE)

n1=10
n2=10
n=n1+n2-2
p=3

y=sqrt((n1*n2)/(n1+n2))*C%*%(x1_bar-x2_bar)

S=((n1-1)*S1+(n2-1)*S2)/(n1+n2+2)

S_new=C%*%S%*%t(C) #pooled S

library(MASS)
T2=t(y)%*%ginv(S_new)%*%y #T2 statistic

F=T2*(n-p+1)/(n*p) #F-statistic
F

qf(.975,p,n-p+1)   #upper 5% point of F

#confidence interval

T.ci <- function(mu, Sigma, n, avec=rep(1,length(mu)), level=0.95){
p <- length(mu)
if(nrow(Sigma)!=p) stop("Need length(mu) == nrow(Sigma).")
if(ncol(Sigma)!=p) stop("Need length(mu) == ncol(Sigma).")
if(length(avec)!=p) stop("Need length(mu) == length(avec).")
if(level <=0 | level >= 1) stop("Need 0 < level < 1.")
cval <- qf(level, p, n-p) * p * (n-1) / (n-p)
zhat <- crossprod(avec, mu)
zvar <- crossprod(avec, Sigma %*% avec) / n
const <- sqrt(cval * zvar)
c(lower = zhat - const, upper = zhat + const)
}

T.ci(mu=C%*%(x1_bar-x2_bar), Sigma=S_new, n=n, avec=c(1,0,0))
T.ci(mu=C%*%(x1_bar-x2_bar), Sigma=S_new, n=n, avec=c(0,1,0))
T.ci(mu=C%*%(x1_bar-x2_bar), Sigma=S_new, n=n, avec=c(0,1,0))

#calculating confidence interval by Bonferroni method

x_bar=C%*%(x1_bar-x2_bar)
TCI <- bon <- NULL
alpha <- 1 - 0.05/(2*3)
for(k in 1:3){
avec <- rep(0, 3)
avec[k] <- 1
TCI <- c(TCI, T.ci(x_bar, S_new, n, avec))
bon <- c(bon,
x_bar[k] - sqrt(S_new[k,k]/n) * qt(alpha, df=n-1),
x_bar[k] + sqrt(S_new[k,k]/n) * qt(alpha, df=n-1))
}
rtab <- rbind(TCI, bon)
round(rtab, 2)