# R-Code zu Blatt 3/Aufgabe 3c)

load("buligoals.rda")

## create differences between seasons
d_i <- goals1314 - goals1415

## number of observations, number of variables
n <- nrow(d_i)
p <- ncol(d_i)


################
### slow (rounded) version

## mean differences per variable
d.bar <- round(colMeans(d_i),1)

## center differences for covariance matrix
centered.d_i <- t(apply(d_i,1,function(x){x-d.bar}))
## estimate covariane matrix
S <- 1/(n-1)*t(centered.d_i)%*%centered.d_i
S <- round(S,1)

## invert covariance matrix
Sinv <- matrix(c(S[4],-S[3],-S[2],S[1]),ncol=2)
factor.Sinv <- (S[1]*S[4]-S[2]*S[3])

## calculate test statistic
T <- ((n-p)*n)/((n-1)*p)/factor.Sinv * (d.bar%*%Sinv%*%d.bar)

## critical region of test
qf(0.95, 2,14)

## p-value
1-pf(T,p,n-p)

################
### fast (exact) version

## mean differences per variable
d.bar2 <- colMeans(d_i)

## estimate covariane matrix
S2 <- cov(d_i)

## or 
cov(goals1314)+cov(goals1415)-cov(goals1314,goals1415)-cov(goals1415,goals1314)

## invert covariance matrix
Sinv2 <- solve(S2)

## calculate test statistic
T2 <- ((n-p)*n)/((n-1)*p) * (d.bar2%*%Sinv2%*%d.bar2)

## critical region of test
qf(0.95, p, n-p)

## p-value
1-pf(T2,p,n-p)