missing}. Note that equation (13.32) reduces to

data matrix, is no greater than min(m, n ’ 1).

± i = x, e i when there are no gaps in x.

Thus SVD extracts the same information from

the sample as a conventional EOF analysis.

13.2.9 Computing Eigenvalues and Eigenvec-

tors. One approach to computing eigenvalues 13.3 Inference

and eigenvectors is to use a ˜canned™ eigen-

analysis routine such as those that are contained 13.3.1 General. We consider the reliability of

in EISPACK [352], IMSL [193], or NAG [298]. eigenvalues and EOF estimates in this section.

Press et al. [322, p. 454] discuss the origins of This is a somewhat different question from that

these routines and give further references. which users generally have in mind when they

enquire about the ˜signi¬cance™ of an EOF. The

An alternative approach uses Singular Value

null hypothesis that is usually implicit in the latter

Decomposition (Appendix B, and see also Press

is that the EOF in question describes only an aspect

et al. [322, pp. 51“63] and Kelly [218, 219]). The

SVD of the conjugate transpose of the m — n of the covariance structure of the ˜noise™ in the

observed system, and the alternative hypothesis is

centred data matrix

that the EOF also describes at least part of the

1

(13.33) dynamics of the observed system. Unfortunately,

X =X I’ J,

n discrimination between ˜noise™ and ˜signal™ in

17 We have implicitly assumed here that m ¤ n. The problem

where X is given by equation (13.29), I is the n—n

is approached similarly when m > n, except we begin by

identity matrix, and J is the n — n matrix of units, obtaining the SVD of X . Note also that in some texts U and

is V are n — n and m — m orthogonal matrices respectively and

S is n — m. The singular values are placed in the diagonal

part of S and the rest of the matrix is zero. We use the

X = USV † ,

†

(13.34)

decomposition given in (13.34) because it is commonly used

in SVD subroutines (see, e.g., Press et al. [322]).

where U is n — m, S and V are each m — m, 18 The rank of a matrix is the dimension of the sub-space

n is the sample size, and m is the dimension of spanned by the columns of that matrix.

13: Empirical Orthogonal Functions

302

13.3.3 The Bias in Estimating Eigenvalues. It

this way is fraught with dif¬culty. We discuss

is natural to ask questions about the reliability

this further in [13.3.4]. However, we ¬rst brie¬‚y

of eigenvalues and EOF estimates, such as the

consider the variance of EOF coef¬cients in

extent to which the estimated patterns resemble

[13.3.2] and the bias of eigenvalue estimates in

the true patterns and how close the estimated and

[13.3.3]. We consider the sampling error of the

true eigenvalues are. These questions do not have

EOFs themselves in [13.3.5,6].

completely satisfactory answers, but there are a

number of potentially useful facts. One of these

13.3.2 The Variance of EOF Coef¬cients of facts is the following set of asymptotic formulae

a Given Set of Estimated EOFs. Assume we that apply to eigenvalue estimates computed from

are given a set of eigenvalues »i and EOFs e i samples that can be represented by n independent

that are derived from a ¬nite sample {x1 . . . xn }. and identically distributed normal random vectors

Then any random vector X can be represented (Lawley [245]):

«

in the space spanned by these estimated EOFs

by using the transformation ± = P † X. The

¬ »j ·

m

¬1 + 1

transformed random variables ± i = X, e i have · + O(n ’2 )

E(»i ) = »i

» ’ »j

their own moments, such as variances Var(± i ) and n

j=1 i

covariances Cov ± i , ± j . In the following we view j=i

the estimated EOFs e i as being ˜¬xed™ (or frozen) (13.38)

«

rather than random.

2 »i ¬ 2·

m

»j

2

¬1 ’ 1

Intuitively one would hope that the variance of ·

Var(»i ) =

n »i ’ » j

± i is equal to that of the real EOF coef¬cient ±i . n

j=1

Unfortunately, this is not the case (see [13.2.6]). j=i

Consider, for example, the ¬rst EOF e 1 and the

+ O(n ’3 ). (13.39)

corresponding EOF coef¬cient ±1 = X, e 1 . The

¬rst EOF minimizes As usual, m is the dimension of the random

vectors. The symbol O(n ’s ) represents a term that

2

1 = E X ’ X, e e .

11

converges to zero as n ’ ∞ at least as quickly as

n ’s does.

By equations (13.38) and (13.39), the eigen-

Replacing e 1 with any other vector, such as e 1 ,

value estimators are consistent:

increases 1 . Thus

lim E (»i ’ »i )2 = 0. (13.40)

Var(X) ’ Var(±1 ) n’∞

2

=E X ’ X, e 1 e 1 However, something unwanted is hidden in

equation (13.38), namely that the estimators of

1 e1 2

< E X ’ X, e the largest and of the smallest eigenvalues are

biased. For the largest eigenvalues, almost all of

= Var(X) ’ Var(± 1 ),

the denominators in equation (13.38) are positive

so that the entire sum is positive, that is, E(»i ) >

that is, Var(±1 ) > Var(± 1 ). Similar arguments lead

»i for large eigenvalues. Similarly, E(»i ) < »i for

to

the smallest eigenvalues.

Together with the results from [13.3.2] this

• Var(± i ) < Var(±i ) for the ¬rst few

EOFs (those corresponding to the largest ¬nding shows that

eigenvalues). • for the largest eigenvalues »i ,

=

Since the total variance Var(X) E(»i ) > »i = Var(±i ) > Var(± i ), (13.41)

m

j=1 Var x j is estimated with nearly zero

bias by Var (X) = m » j , it follows that • for the smallest eigenvalues »i ,

j=1

• Var(± i ) > Var(±i ) for the last few EOFs. E(»i ) < »i = Var(±i ) < Var(± i ). (13.42)

Examples show that these deviations may be Relations (13.41) and (13.42) illustrate that

considerable, in particular for small eigenvalues we must be cautious when using estimated

[392]. eigenvalues. First, the estimates are biased: the

13.3: Inference 303

large eigenvalues are overestimated and the small and noise components (recall [10.1.1]), has

eigenvalues »1 , . . . , »m and EOFs e 1 , . . . , e m .

ones are underestimated. More important, though,

is the inequality E(»i ) = Var(± i ): The sample Now construct a multivariate white noise Zt from

eigenvalue » is a biased estimator of the variance iid N (0, I) random vectors. Then the multivariate

of ± i = white noise process Yt = P 1/2 Zt , which is

X, e i , for any frozen set of

estimated EOFs e i . Similarly, Cov ± i , ± j = completely devoid of ˜dynamics,™ has the same

Cov(± i , ± j ) = 0. eigenvalues and eigenvectors as Xt . Thus we

cannot always diagnose dynamical structure from

the zero lag covariance structure of a process.21

13.3.4 ˜Selection Rules.™ Many so-called selec-

Our recommendation is to avoid using selection