MAST90138 Week 3 Lab
• Recall that if Z1, . . . , Zn are independent N(0, 1) then
n
X = Z k2 ∼ χ 2n
k=1
is a chi square with n degrees of freedom.
• Recall that if M is an p × n matrix whose columns are independent and all have a Np(0, Σ) distri-
bution, then
Problems
Y=MMT ∼Wp(Σ,n). (1)
1. In R, generate a sample X1, . . . , Xn of size n = 200 from a N2(μ, Σ) distribution and draw the scatterplot of the pairs (X1i, X2i), for i = 1, . . . , n, using a red ∗ as the symbol to represent the data points. Do this for the following four different matrices Σ:
10 51 52 52 Σ=02,Σ=12,Σ=22,Σ=21
and plot the four scatterplots on a single graph window showing 4 scatterplots presented as the first two scatterplots on the first row and the last two scatterplots on the second row. Label the axes as X1andX2andaddthetitleScatterplot 1,…,Scatterplot 4tothefirst,…,fourthscatterplot. What do you notice about the cloud of points when you compare the four cases? Explain what is going on.
2. Show that if X ∼ Np(μ, Σ) and Σ is invertible, then
Y =(X−μ)TΣ−1(X−μ)∼χ2p. (2)
3. Let p = 1 and take Σ = σ2, a number. Show that Y defined above is W1(σ2, n), is also equal to σ2 times a χ2n.
4. Show that if Y is defined at (1) and B is a q×p matrix then BYBT ∼Wq(BΣBT,n).
5. ShowthatifY isdefinedat(1)andaisap×1vectorsuchthataTΣa̸=0,then aT Ya/aT Σa ∼ χ2n .
1
Solution to 1:
We see that in the first graph, the data seems to create some sort of a disk. This is because the two components of the vector are independent (as their covariance is equal to zero) and so the data are not located around an oblique line. Recall that the contour lines of a bivariate normal distribution are ellipses, and in the case where the components are independent, they are circles, whence the shape of the scatterplot. In the other three examples, the two components of the normal vector are not independent as the contour lines of their distribution become ellipses. From case 2 to 4, the ellipses become increasingly more tight around an oblique line because the correlation between the two components increases. Recall that we can compute the correlation matrices by taking:
R = diag(σ11, σ22)−1/2Σ diag(σ11, σ22)−1/2.
We deduce that from the first to the fourth case, ρ12 = 0, 0.3162278, 0.6324555, 0.8944272
library(MASS)
mu = c(1, 2)
sigma1 = matrix(c(1, 0, 0, 2), nr= 2, byrow = T)
sigma2 = matrix(c(5, 1, 1, 2), nr= 2, byrow = T)
sigma3 = matrix(c(5, 2, 2, 2), nr= 2, byrow = T)
sigma4 = matrix(c(5, 2, 2, 1), nr= 2, byrow = T)
N = 200
X1 = mvrnorm(N, mu = mu , Sigma = sigma1)
X2 = mvrnorm(N, mu = mu , Sigma = sigma2)
X3 = mvrnorm(N, mu = mu , Sigma = sigma3)
X4 = mvrnorm(N, mu = mu , Sigma = sigma4)
par(mfcol = c(2, 2))
plot(X1[,1], X1[,2], xlab = “X1”, ylab = “X2″, main =”Scatterplot 1″, pch =”*”, col =”red”)
plot(X2[,1], X2[,2], xlab = “X1”, ylab = “X2″, main =”Scatterplot 2″, pch =”*”, col =”red”)
plot(X3[,1], X3[,2], xlab = “X1”, ylab = “X2″, main =”Scatterplot 3″, pch =”*”, col =”red”)
plot(X4[,1], X4[,2], xlab = “X1”, ylab = “X2″, main =”Scatterplot 4″, pch =”*”, col =”red”)
D1 = diag(diag(1/sqrt(sigma1)))
D2 = diag(diag(1/sqrt(sigma2)))
D3 = diag(diag(1/sqrt(sigma3)))
D4 = diag(diag(1/sqrt(sigma4)))
R1 = D1%*%sigma1 %*%D1
R2 = D2%*%sigma2 %*%D2
R3 = D3%*%sigma3 %*%D3
R4 = D4%*%sigma4 %*%D4
2