## October 25, multivariate Q-Q plots for assesment of joint normality.
## Also a full re-visit of the glucose data
rm(list = ls())

###################################################################################
## example from the book
y <- read.table("glucose.txt")
n <- nrow(y)
y <- as.matrix(y)
y[1:3,]

## Univariate
############################
## Boxplots
boxplot(y, col="yellow")

## normal QQ plots
par(mfrow=c(2,3))
for(i in 1:6){
qqnorm(y[,i], pch=16)
qqline(y[,i])
} 


## Bivariate
pairs(y)

## Multivariate
Mu.est <- apply(y, 2, mean)
Sigma.est <- cov(y)

d.sq <- rep(0, n)
for(i in 1:n){
d.sq[i] <- (y[i,] - Mu.est)%*%solve(Sigma.est)%*%(y[i,] - Mu.est)  ## make sure to becareful about dimensions.  This is correct here.
}
d.sq.sort <- sort(d.sq)


# get the quantiles of the Chi-Squared distribution with p=6 degrees of freedom
par(mfrow=c(1,1))
qu.chi <- qchisq( (1:n-0.5)/n, df=6)
plot(qu.chi, d.sq.sort,pch=16, xlab="Theoretical Quantiles", ylab="Sample Quantiles", main="Chi Squared Q-Q plot")

## let's add a line which passes throuugh the first and third quartile
q.25 <- qchisq(0.25, df=6)
q.75 <- qchisq(0.75, df=6)
d.sq.25 <- quantile(d.sq, prob=0.25)
d.sq.75 <- quantile(d.sq, prob=0.75)

b <- (d.sq.75-d.sq.25)/(q.75-q.25)
a <- d.sq.25 - b*q.25
abline(a,b,lwd=2)