2
See Exercise 5.17 to show that the Tp,m distribution does not depend on the choice
2
of ⌃. If we take p = 1, Hotelling’s T distribution is equal to the square of a student
tm distribution. What if p > 1? Well, then it is linked to the F distribution, which
we now recall:
Definition 5.1.23 (F distribution). Let p, q 2 {1, 2, . . .}. A random
variable W 2 R is said to follow an Fp,q distribution, also called an F
distribution with (p, q) degrees of freedom, if
U/p
W = ,
V /q
2 2
for some independent U ⇠ p, V ⇠ q.
d
Notice in particular that (tq )2 = F1,q , therefore the F distribution generalises the stu-
dent’s t distribution. It turns out that a re-scaled version of a Hotelling’s T 2 follows
an F distribution. Specifically, we have the following result.
Proposition 5.1.24 (Link between Hotelling’s and the F distribu-
2
tion). If D ⇠ Tp,m , then
m p+1
D ⇠ Fp,m p+1 .
mp
Hence, we can easily evaluate the density and tail probabilities using those from the
F distribution. The proof of this Proposition is beyond the scope of this module.
5.2 Parameter estimation
We now discuss how to obtain maximum likelihood estimates for the mean and
covariance of a multivariate Normal distribution, and their associated confidence
regions. We will only consider the case where ⌃ is invertible. If it is not, then
one needs to work in coordinates with respect to the eigenvectors of ⌃ that have
non-zero eigenvalues (and this is beyond the scope of the module).
5.2.1 Point estimates
Let x1 , . . . , xn be i.i.d. draws from a Np (µ, ⌃). The likelihood function is
n ⇢ !
Y 1 1
L(µ, ⌃) = p 1 exp (xi µ)T ⌃ 1 (xi µ)
i=1 (2⇡) det(⌃)
2 2 2
( n
)
1 1X
= np n exp (xi µ)T ⌃ 1 (xi µ) (5.2.1)
(2⇡) 2 det(⌃) 2 2 i=1
The maximum likelihood estimate (MLE) is obtained by maximising (5.2.1) with
respect to (µ, ⌃). In order to find the MLE we will need the following result.
96
Lemma 5.2.1 (MLE technical lemma). Let B be a symmetric positive
definite (SPD) p ⇥ p matrix and b > 0. Then the maximum of the function
✓ ◆
1
g(⌃) , det(⌃) exp
b 1
Tr(⌃ B) ,
2
1
defined over the space of SPD matrices, is obtained for ⌃ = 2b
B.
Proof. Non-examinable.
Proposition 5.2.2 (MLE for the mean and covariance of MVN dis-
tribution). Let x1 , . . . , xn be an observed random sample from Np (µ, ⌃).
1. The MLE for µ is
n
1X
µ̂ , x̄ , xi
n i=1
2. The MLE for ⌃ is
X n
ˆ, 1
⌃ (xi x̄)(xi x̄)T .
n i=1
Proof. We start by maximising (5.2.1) with respect to µ, which is equivalent to
minimising
Xn
(xi µ)T ⌃ 1 (xi µ).
i=1
The optimum could be easily found using vector derivatives, but here we shall use
a more basic approach. Using that the trace of a number is that same number and
other properties of the trace, we re-write the expression as
n
! n
X X
Tr (xi µ)T ⌃ 1 (xi µ) = Tr (xi µ)T ⌃ 1 (xi µ)
i=1 i=1
n n
!
X X
= Tr ⌃ 1 (xi µ)(xi µ)T = Tr ⌃ 1
(xi µ)(xi µ)T
i=1 i=1
n
!
X
= Tr ⌃ 1
(xi x̄ + x̄ µ)(xi x̄ + x̄ µ)T
i=1
n
!
X T
1
= Tr ⌃ (xi x̄) + (x̄ µ) (xi x̄) + (x̄ µ)
i=1
n n
!
X X
1 T 1 T
= Tr ⌃ (xi x̄)(xi x̄) + ⌃ (x̄ µ)(x̄ µ) (5.2.2)
i=1 i=1
n
!
X
= Tr ⌃ 1
(xi x̄)(xi x̄)T + n(x̄ µ)T ⌃ 1 (x̄ µ). (5.2.3)
i=1
97
Since the first term in (5.2.3) does not depend on µ, we just need to minimise
(x̄ µ)T ⌃ 1 (x̄ µ). This expression is greater than or equal to zero given that ⌃ 1
is positive definite, and is exactly equal to zero at µ = x̄. Therefore the minimum
is achieved at µ = x̄.
To find the MLE for ⌃, we plug in µ = x̄ as this maximises the likelihood for
any ⌃. Taking (5.2.1) and (5.2.2) the goal is to maximise
1 1
Tr(⌃ 1B
) = g(⌃),
e 2
det(⌃)n/2
Pn
where g(·) is defined in Lemma 5.2.1, with b = n/2 and B = i=1 (xi x̄)(xi x̄)T .
Lemma 5.2.1 tells us that the maximum occurs at
n
1 1X
nB = (xi x̄)(xi x̄)T .
22 n i=1
Similarly to the univariate Normal case, the MLE for ⌃ is biased and typically
one uses the unbiased version with n 1 in the denominator, that is
n
X
1
S= (Xi X)(Xi X)T
n 1 i=1
(recall Section 4.2.4.)
5.2.2 Confidence regions
1
Pn
The following result gives the sampling distribution for X and S = n 1 i=1 (Xi
X)(Xi X)T for Normally distributed data.
Proposition 5.2.3. Let X1 , . . . , Xn be an i.i.d. sample from Np (µ, ⌃).
P
1. X = n1 ni=1 Xi ⇠ Np (µ, n1 ⌃).
2. (n 1)S ⇠ Wp (⌃, n 1).
3. X and S are independent.
Proof. 1. Here is a sketch, the details of which you should complete on your own.
Since the Xi s are i.i.d. MVN, they are jointly MVN, and therefore the sample
mean (being a linear transformation of a multivariate normal) is MVN. One di-
rectly computes
⇥ the expectation⇤ of X. For its covariance, we use the formula
Cov(Y) = E Y(Y E Y) (see Exercise 4.21) and use the linearity of the ex-
T
pectation.
2. Let X be the n ⇥ p matrix with XTi as its i-th row. Recall that if H =
(I 11T /n), then HX returns the column-centered data matrix, and hence since
98
HH = H we have S = XT HX/(n 1). Using Q2)c) from Assignment 1 we therefore
have
n 1
X
(n 1)S = (XT uj )(XT uj )T .
j=1
for some orthonormal u1 , . . . , un 1 with uTj 1 = 0 for all j. Now we know that the
vectors (XT uj )j are jointly multivariate normal, because they are linear transforma-
tions of X, which is MVN. Some calculations (similar to those in Assignment 1) give
that that they are i.i.d. Np (0, ⌃). We then just apply the definition of the Wishart
distribution to conclude the proof of part 2.
3. Notice that (n 1)S = (HX)T (HX), so it is enough to show that HX is inde-
pendent of X. Since they are jointly Gaussian (convince yourself of this), Proposi-
tion 5.1.9 tells us that it is sufficient to check that HX and X are uncorrelated. We
can do this by proving that X is uncorrelated with each row of HX, given by Xi X
for i = 1, . . . , n. But we do indeed have
Cov(X, Xi X) = E{X(Xi X)T } = E{X(Xi µ)T } E{X(X µ)T }
1 1 1
= E{Xi (Xi µ)T } Cov(X) = ⌃ ⌃ = 0,
n n n
where we have used Exercise 4.21 and part 1 of this result.
Example 5.2.4. Let X = (X1 , X2 , X3 )T ⇠ Np (µ, ⌃) with µ = (10, 5, 0)T and
0 1
2 1 1
⌃ = @1 2 1 A .
1 1 2
Define Y1 = X1 X2 and Y2 = X3 and suppose we have a sample of size n = 20
from X. What is the sampling distribution of Ȳ?
We first note that Y = (Y1 , Y2 )T = CX, with
✓ ◆
1 1 0
C= .
0 0 1
Therefore the distribution of✓Y is bivariate
◆ Normal with mean µy = Cµ = (5, 0)T
2 0
and covariance ⌃y = C⌃C T = . Hence, the distribution of the sample mean
0 2
based on n = 20 observations is
Ȳ ⇠ N2 (µy , ⌃y /20).
The previous result tells us the sampling distribution of the MLE (X, nn 1 S) for
Normally distributed data.
In Proposition 5.1.18 we saw that, for X ⇠ Np (µ, n1 ⌃), the statistic
n(X µ)T ⌃ 1 (X µ)
99
follows a 2p distribution. Hence, if ⌃ were known we could use this result to obtain
a 1 ↵ confidence region for µ as follows:
µ 2 Rp : n(X µ)T ⌃ 1 (X µ) 2
p (1 ↵)
When ⌃ is unknown, we can replace it with an estimate and use the following result.
Proposition 5.2.5 (Confidence regions for µ).
iid P
Let X1 , . . . , Xn ⇠ Np (µ, ⌃), and let X = n1 ni=1 Xi and
n
X
1
S= (Xi X)(Xi X)T .
n 1 i=1
If ⌃ is invertible, and n p + 1, then
(n 1)p
n(X µ)T S 1
(X µ) ⇠ Fp,n p
n p
In particular, a 1 ↵ confidence region for µ is given by the ellipsoid
⇢
(n 1)p
y 2 Rp : n(X y)T S 1 (X y) Fp,n p (1 ↵) ,
n p
where Fp,n p (1 ↵) is the 1 ↵ quantile of an Fp,n p distribution.
p p
Proof. Define T 2 = n(X µ)T S 1 n(X µ), and note that this is the product
of a Np (0, ⌃) times an independent ( n 1 1 Wp (⌃, n 1)) 1 times the same Np (0, ⌃).
By Definition 5.1.22, T 2 follows Hotelling’s distribution with degrees of freedom
(p, n 1), which gives the result after applying Proposition 5.1.24.
A basic result in probability theory is that when n ! 1, (nn 1)p
p
Fp,n p converges
2
to a p . Intuitively this makes sense, as we know by the law of large numbers
p
and Slutsky’s Lemma that S ! ⌃ as n ! 1, and thus applying the continuous
d
mapping theorem we get T 2 ! 2p .
Example 5.2.6. We continue Example 5.1.19, but just use 20 draws from a bivariate
normal. We compute the sample mean X = (0.3, 0.06) and sample covariance
✓ ◆
1.51 0.97
S= .
0.97 1
We compute 95% confidence regions for µ, both using the true ⌃ and its estimate
S. When ⌃ is known, the confidence region is given by the ellipse
{y 2 R2 : (X y)T ⌃ 1 (X y) 2
2 (0.95)/20}.
When ⌃ is unknown, we are looking for the ellipse
Fp,18 (0.95) 19 ⇥ 2
{y 2 Rp : (X y)T S 1
(X y) ⇥ }
20 18
100
These ellipses can be easily plotted in R using the ellipse function and indicat-
ing the cuto↵ point with the argument t. The R code is below and Figure 36 shows
the results. On the upper panel we see the original draws and that both confidence
regions are quite similar. On the lower panel we provide a zoom-in where we can
appreciate that the 95% confidence region using S is slighly larger than that for ⌃
(as it’s usually, though not always, the case). Both regions are centered at X and
contain the true value µ = (0, 0)T .
The R code below shows that the critical values for the known and unknown cases
are not quite similar (0.3 vs. 0.375), as expected since n is not large. In particular,
this means that the confidence region without knowing the true covariance is larger
to the confidence region when the true covariance is known. This di↵erence becomes
smaller if n is larger.
op <- par(mfrow=c(1,1), mai=rep(.4,4))
plot(x,xlab='X1',ylab='X2', xlim=c(-4,4), ylim=c(-2,2))
lines(ellipse(Sest, centre=muest, t=sqrt(tf)), lwd=2)
lines(ellipse(S, centre=muest, t=sqrt(tchi)),col=2,lty=2, lwd=2)
points(muest[1],muest[2], pch=15, cex=1.5)
points(mu[1],mu[2], pch=17, cex=1.5, col=4)
#text(muest[1],muest[2],'Sample mean',pos=1,cex=1)
#text(mu[1],mu[2],'True mean',pos=1,cex=1,col='blue')
legend('topright', legend=c('True mean', 'Sample mean'),
pch=c(17, 15), col=c(4,1))
par(op)
Proposition 5.2.5 tells us how to obtain confidence regions for the mean of a single
population. Often we will have data arising from several populations, in which case
we may be interested in obtaining confidence regions for each group. In general, we
can simply apply Proposition 5.2.5 to each group separately. However, when the
covariances in each group are equal (or in practice, similar enough) and the sample
size per group is limited, it may be preferrable to obtain an overall estimate for ⌃
across all groups. In that case the following result applies.
Proposition 5.2.7 (Confidence regions for µ1 , . . . , µK ).
Let X1 , . . . , XK be the sample means from K independent multivariate
Normal populations with means µ1 , . . . , µK , and common covariance ⌃.
Let
XK
1
Sp = (ni 1)Si ,
N K i=1
where n1 , . . . , nK and S1 , . . . , SK are the sample sises
PK and sample covariance
matrices in each group (respectively) and N = i=1 ni . Then
(N K)p
ni (Xi µi )T Sp 1 (Xi µi ) ⇠ Fp,N K p+1
N K p+1
101
True mean
2
Sample mean
●
●
● ●
1
● ●
●
●
●
●
● ●
0
●
●
●
●
−1
●
−2
●
−4 −2 0 2 4
Figure 36: Bivariate Normal data and 95% confidence region using true (dashed
red) and estimated covariance (solid black).
102
In particular, a 1 ↵ confidence region for µi is given by the ellipse
⇢
(N K)p
yi 2 Rp : ni (Xi yi )T Sp 1 (Xi yi ) Fp,N K p+1 (1 ↵)
N K p+1
Proof. Because of the additive property of independent Wishart distributions, the
pooled sample covariance Sp follows a Wishart with N K degrees of freedom.
Hence the result follows immediately from the proof of Proposition 5.2.5, by simply
noting that we have the product of a Np (0, ⌃) times the inverse of an independent
1
N K
Wp (⌃, N K) times the same Np (0, ⌃).
5.2.3 Asymptotics of the Sample Mean
We have see that X is the MLE of the mean of a MVN random vector. However
(provided ⌃ is finite), X may be a sensible estimator for non-normally distributed
data as well, as it estimates the population mean consistently as n ! 1. Also,
when n ! 1 then X is approximately Normal even when the individual Xi ’s are
not Normally distributed.
Proposition 5.2.8 (Law of Large Numbers and Central Limit The-
orem).
iid
Let X1 , . . . , Xn ⇠ X be an i.i.d. sample from a multivariate distribution.
1. If E|X| < 1, then
a.s.
X ! µ, as n ! 1,
where µ = E X.
2. If E|X|2 < 1, then
p d
n(X µ) ! Np (0, ⌃), as n ! 1,
where ⌃ = Cov(X).
Proof. Part (1) follows from applying the univariate Strong Law of Large Numbers
for each element in X.
We do not prove Part (2) here, but the result follows from seeing that the char-
acteristic function of X converges to that of a multivariate Normal (analogously to
the proof for the univariate Central Limit Theorem).
5.3 Hypothesis testing
We now focus on testing hypotheses, that is on determining whether a given param-
eter value (or function of parameter values) is likely to have generated the observed
data. These multivariate tests generalise the usual one-sample t-test, two-sample
t-test and F-test (ANOVA) to the case where there are several response variables.
103