Multivariate Analysis

Agricultural Statistics & Data Science · Statistical theory core · also Multivariate Statistical Methods, Applied Multivariate Analysis

12units
12chapters
5,267words
102practice questions

Analysing many correlated variables together: the multivariate normal and its estimation, Wishart and Hotelling's T-squared, MANOVA, and the reduction and classification tools agricultural scientists use daily (principal components, factor analysis, discriminant and cluster analysis, Mahalanobis D-squared).

In the exams: AICE-JRF Ag Stat Unit 2 (heavy: MVN, Wishart, Hotelling, MANOVA, PCA, FA, discriminant, cluster named); GATE ST 'Multivariate Analysis' (heavy); UPSC ISS Statistics II (heavy); GATE DA PCA (present)

Syllabus

  1. Multivariate normal distribution: properties, marginals and conditionals
  2. Estimation of the mean vector and dispersion matrix
  3. Wishart distribution
  4. Hotelling's T-squared and its applications
  5. Multivariate analysis of variance (MANOVA)
  6. Principal component analysis
  7. Factor analysis
  8. Canonical correlation analysis
  9. Discriminant analysis and classification
  10. Cluster analysis: hierarchical and k-means
  11. Mahalanobis D-squared and genetic divergence
  12. Multidimensional scaling and correspondence analysis

Study manual

Multivariate analysis studies several correlated variables measured on the same units, treating them jointly rather than one at a time. A plant breeder records yield, days to flowering, plant height, grain weight and protein on every line; a soil scientist measures a dozen nutrients on every sample; a survey records many characteristics of every farm. Analysing each variable separately ignores their correlations, multiplies the chance of false findings and misses patterns that exist only in combination. Multivariate methods test hypotheses about whole vectors of means, reduce many variables to a few informative combinations, and classify or group units on all variables at once.

The subject has two halves. The inferential half rests on the multivariate normal distribution: estimation of the mean vector and dispersion matrix, the Wishart distribution of the sample sum-of-squares-and-products matrix, Hotelling's T² for mean vectors and multivariate analysis of variance (MANOVA). The exploratory half contains the tools agricultural scientists use daily: principal components, factor analysis, canonical correlation, discriminant analysis, cluster analysis, Mahalanobis D² for genetic divergence, multidimensional scaling and correspondence analysis.

The manual builds on Linear Models and Regression (matrix algebra of least squares, multiple and partial correlation), Probability Distributions and Sampling Distributions (the bivariate normal, chi-square and F, quadratic forms and Cochran's theorem) and Statistical Inference I and II (maximum likelihood and likelihood-ratio tests), and refers to them instead of repeating them. Matrix notation is used throughout: vectors are columns, A′ is the transpose, |A| the determinant, tr(A) the trace and Σ the population dispersion (covariance) matrix.

Examinations test it heavily. The ICAR AICE-JRF/SRF Agricultural Statistics paper names the multivariate normal, Wishart, Hotelling's T², MANOVA, principal components, factor, discriminant and cluster analysis in Unit 2; the Genetics and Plant Breeding papers ask about D² and genetic divergence. GATE Statistics has a multivariate section, GATE Data Science asks about PCA, and the Indian Statistical Service covers the theory in Statistics II.

1. The multivariate normal distribution

Definition and basic properties

A p-dimensional random vector X has a multivariate normal distribution Nₚ(μ, Σ), with mean vector μ and positive definite dispersion matrix Σ, if its density is f(x) = (2π)^(−p/2)|Σ|^(−1/2) exp[−(1/2)(x − μ)′Σ⁻¹(x − μ)]. An equivalent definition, which also covers singular Σ, is that every linear combination a′X is univariate normal. The bivariate normal of the distributions manual is the case p = 2.

Key properties. Linear transformations stay normal: AX + b ~ N(Aμ + b, AΣA′). Every subset of the variables is multivariate normal, with the corresponding sub-vector of μ and sub-matrix of Σ. Two sub-vectors are independent if and only if their cross-covariance matrix is zero, the property that makes zero correlation equivalent to independence under joint normality. The quadratic form (X − μ)′Σ⁻¹(X − μ) is χ²ₚ, so the contours of constant density are ellipsoids, and the ellipsoid (x − μ)′Σ⁻¹(x − μ) = χ²ₚ,α contains probability 1 − α. Normal marginals do not imply joint normality.

Conditional distributions

Partition X into X₁ (q components) and X₂ (p − q components), with μ and Σ partitioned conformably into μ₁, μ₂, Σ₁₁, Σ₁₂, Σ₂₁ and Σ₂₂. Then the conditional distribution of X₁ given X₂ = x₂ is normal with mean μ₁ + Σ₁₂Σ₂₂⁻¹(x₂ − μ₂) and dispersion Σ₁₁ − Σ₁₂Σ₂₂⁻¹Σ₂₁. The conditional mean is linear in x₂, which is why regression on normally distributed predictors is linear, and the conditional dispersion does not depend on x₂ (homoscedasticity). The matrix Σ₁₂Σ₂₂⁻¹ contains the regression coefficients, and Σ₁₁·₂ = Σ₁₁ − Σ₁₂Σ₂₂⁻¹Σ₂₁ is the partial dispersion matrix from which partial correlations are computed.

Example: yield (X₁) and two predictors, plant height and tiller number (X₂), are jointly normal. The best predictor of yield from height and tillers is the conditional mean, a linear function with coefficients Σ₁₂Σ₂₂⁻¹, and the proportion of yield variance explained is the squared multiple correlation Σ₁₂Σ₂₂⁻¹Σ₂₁/σ₁₁. Selection indices in plant and animal breeding (the Smith-Hazel index) use exactly this algebra: the index weights are b = P⁻¹Ga, with P the phenotypic and G the genetic dispersion matrices.

Checking multivariate normality in practice uses univariate normal plots of each variable, scatter plots of pairs, and a chi-square plot of the squared Mahalanobis distances of observations from the mean, which should lie close to a straight line against χ²ₚ quantiles. Mardia's measures of multivariate skewness and kurtosis provide formal tests.

Remember

2. Estimating the mean vector and dispersion matrix

Sample estimates and their distributions

From a random sample x₁, ..., xₙ of p-variate observations, the sample mean vector is x̄ = (1/n)Σxᵢ and the sample dispersion matrix is S = [1/(n − 1)]Σ(xᵢ − x̄)(xᵢ − x̄)′. The matrix A = (n − 1)S is the corrected sum-of-squares-and-products (SSCP) matrix: its diagonal holds the corrected sums of squares of each variable and its off-diagonal the corrected sums of products. The sample correlation matrix R has entries rⱼₖ = sⱼₖ/√(sⱼⱼsₖₖ).

Under multivariate normality: x̄ ~ Nₚ(μ, Σ/n); A = (n − 1)S has the Wishart distribution Wₚ(n − 1, Σ) (Chapter 3); and x̄ and S are independent, extending the univariate result. S is unbiased for Σ. The maximum-likelihood estimators are μ̂ = x̄ and Σ̂ = A/n = (n − 1)S/n, biased like the univariate divisor-n variance (see the estimation manual for why the MLE of a variance is biased). By invariance, the MLE of the correlation matrix is R. (x̄, S) is sufficient and complete for (μ, Σ).

For S to be non-singular (invertible) the sample size must exceed the number of variables: n > p, and in practice n should be several times p for stable estimates. When p is large relative to n, as in genomic and spectral data, S is singular or nearly so, and regularised (shrinkage) estimators of Σ are used.

Generalised variance and practical points

Two scalar summaries of dispersion are common. The total variance tr(S) is the sum of the variances; the generalised variance |S| is the determinant, proportional to the squared volume of the scatter ellipsoid. |S| is small when variables are highly correlated (the ellipsoid is flat) even if each variance is large, so it summarises both spread and redundancy.

Missing values, outliers and scale differences cause most practical problems. Pairwise deletion can give an S that is not positive definite; the EM algorithm or multiple imputation is preferred. Robust estimators (minimum covariance determinant) resist outliers, which can distort every multivariate method. Variables measured on very different scales (yield in kg/ha, plant height in cm, protein in %) dominate covariance-based methods according to their units, which is why analyses are often based on the correlation matrix, i.e. on standardised variables.

A small worked case shows the arithmetic. Two traits measured on 11 plants give corrected sums of squares 40 and 90 and a corrected sum of products 30. Then S has variances 4 and 9 and covariance 3, the correlation is 3/√36 = 0.5, tr(S) = 13 and |S| = 4 × 9 − 3² = 27. If the correlation had been 0.9 (covariance 5.4), the total variance would be unchanged at 13 but |S| would fall to 36 − 29.16 = 6.84, showing how the generalised variance shrinks as the variables become redundant. The MLE of Σ divides the SSCP matrix by 11 instead of 10.

Remember

3. The Wishart distribution

Definition and properties

If Z₁, ..., Z_m are independent Nₚ(0, Σ) vectors, the p × p random matrix W = ΣZᵢZᵢ′ has the Wishart distribution Wₚ(m, Σ) with m degrees of freedom. It is the multivariate generalisation of the chi-square: for p = 1, W = Σzᵢ² ~ σ²χ²_m. Its mean is E(W) = mΣ. When m ≥ p and Σ is positive definite, W is positive definite with probability one and has a density proportional to |W|^((m − p − 1)/2) exp[−tr(Σ⁻¹W)/2].

Properties used constantly: independent Wishart matrices with the same Σ add, with degrees of freedom adding; for a fixed q × p matrix C, CWC′ ~ W_q(m, CΣC′); in particular, for a fixed vector a, a′Wa/a′Σa ~ χ²_m, so every linear combination of the variables has a chi-square sum of squares; diagonal elements satisfy wⱼⱼ/σⱼⱼ ~ χ²_m. The sample SSCP matrix from n normal observations is Wₚ(n − 1, Σ), one degree of freedom being lost to the mean, exactly as in the univariate case.

The generalised variance has a known distribution: |W|/|Σ| is distributed as the product of p independent chi-squares with m, m − 1, ..., m − p + 1 degrees of freedom (from Bartlett's decomposition of W into a triangular matrix with independent chi and normal entries). Wishart matrices appear in the null distributions of Hotelling's T², Wilks' Λ and the other MANOVA statistics, and the inverse Wishart is the conjugate prior for Σ in Bayesian multivariate analysis.

Using the Wishart in practice

Two practical consequences follow. Because the SSCP matrices for treatments and error in a multivariate analysis of variance are independent Wishart matrices under normality, their ratios behave like multivariate F-ratios, which is how MANOVA test statistics are built (Chapter 5). And because a′Wa/a′Σa is χ², confidence intervals and tests for the variance of any chosen index (a selection index, a ratio of nutrient totals) follow from univariate chi-square theory. The Wishart is sensitive to non-normality in the same way that the chi-square test for a variance is, so inferences about dispersion matrices (Box's M test for equal covariance matrices, tests of sphericity) are much less robust than inferences about means.

In genetics, the Wishart framework underlies estimation of genetic and phenotypic dispersion matrices from multivariate analyses of variance and covariance of breeding trials: the mean-product matrices for genotypes and error are independent Wishart matrices, and their difference, divided by the number of replicates, estimates the genetic dispersion matrix. Because such a difference need not be positive definite in small trials, estimated genetic correlations can exceed 1 in absolute value, a known artefact that signals imprecise estimation rather than a real relationship; restricted maximum likelihood with constraints, or Bayesian estimation with an inverse-Wishart prior, avoids it.

Remember

4. Hotelling's T²

One-sample T²

To test H₀: μ = μ₀ for a p-variate normal mean with unknown Σ, Hotelling's statistic is T² = n(x̄ − μ₀)′S⁻¹(x̄ − μ₀), the multivariate analogue of the squared one-sample t. Under H₀, [(n − p)/(p(n − 1))]T² ~ F(p, n − p). For p = 1 it reduces to t², with F(1, n − 1). T² is the likelihood-ratio test, invariant under non-singular linear transformations of the variables, and it is the UMP invariant test. It also equals the maximum over all linear combinations a of the squared univariate t for a′x: the test asks whether any combination of the variables departs from its hypothesised mean, which is why simultaneous (Scheffé-type) intervals for all a′μ follow from it.

Example: five nutrient concentrations are measured on 20 soil samples from a field, and H₀ is that the mean vector equals a recommended standard. If T² = 18.5, then F = (20 − 5)/(5 × 19) × 18.5 = 2.92 with (5, 15) df, which exceeds F₀.₀₅(5, 15) = 2.90 marginally: the field's nutrient profile differs from the standard, although no single nutrient need differ significantly on its own.

Two-sample T² and applications

For two independent samples from normal populations with a common Σ, the pooled dispersion matrix is S_p = [(n₁ − 1)S₁ + (n₂ − 1)S₂]/(n₁ + n₂ − 2), and T² = [n₁n₂/(n₁ + n₂)](x̄₁ − x̄₂)′S_p⁻¹(x̄₁ − x̄₂). Under H₀: μ₁ = μ₂, [(n₁ + n₂ − p − 1)/(p(n₁ + n₂ − 2))]T² ~ F(p, n₁ + n₂ − p − 1). The quantity (x̄₁ − x̄₂)′S_p⁻¹(x̄₁ − x̄₂) is the sample Mahalanobis D² between the two groups (Chapter 11), so T² is D² scaled by the effective sample size.

Applications include comparing two varieties on several quality traits at once, comparing soil profiles under two managements, paired multivariate comparisons (T² on the vector of differences, with n − 1 error df), repeated measures over time on the same plants (profile analysis tests whether mean profiles are parallel, coincident and flat, using T² on contrasts), and confidence regions for mean vectors. When covariance matrices differ, Behrens-Fisher-type corrections (Yao, Nel and van der Merwe) replace the pooled test. If the assumption of normality is doubtful, T² is fairly robust to non-normality for moderate samples with equal sizes, but sensitive to unequal covariance matrices when sample sizes differ.

Remember

5. Multivariate analysis of variance

SSCP partition and Wilks' lambda

MANOVA extends ANOVA to vector responses. For a one-way layout with g groups and N observations in total, the total SSCP matrix splits into a between-groups matrix B (hypothesis, g − 1 df) and a within-groups matrix W (error, N − g df): T = B + W, the matrix version of the ANOVA identity. Under normality with a common Σ, W ~ Wₚ(N − g, Σ), and under H₀ (equal mean vectors) B ~ Wₚ(g − 1, Σ), independent of W.

Wilks' lambda, the likelihood-ratio statistic, is Λ = |W|/|B + W|, between 0 and 1; small values reject H₀. It equals Π 1/(1 + λᵢ), where λᵢ are the eigenvalues of W⁻¹B. Bartlett's approximation, −[N − 1 − (p + g)/2] ln Λ ≈ χ²_(p(g−1)), is used for large samples; exact F transformations exist for p = 1 or 2 and for g = 2 or 3. For g = 2, Λ is a one-to-one function of the two-sample T². Example: with p = 3 traits, g = 4 varieties and N = 40 plots, a computed Λ = 0.45 gives −[39 − 3.5] ln 0.45 = 35.5 × 0.799 = 28.3, against χ²₉,₀.₀₅ = 16.92: the varieties differ in their mean vectors.

Other test statistics and follow-up

Three other statistics are functions of the eigenvalues λᵢ of W⁻¹B: Pillai's trace V = tr[B(B + W)⁻¹] = Σλᵢ/(1 + λᵢ); the Hotelling-Lawley trace U = tr(W⁻¹B) = Σλᵢ; and Roy's largest root, based on the largest λᵢ. When there is only one non-zero eigenvalue (g = 2 or p = 1) all four are equivalent. Otherwise they differ in power: Roy's root is most powerful when group means differ along a single direction, Pillai's trace is the most robust to violations of assumptions and is often recommended for routine use, and Wilks' Λ is intermediate.

Assumptions are multivariate normality, independence and equal covariance matrices across groups, the last checked with Box's M test (which is itself very sensitive to non-normality). After a significant MANOVA, follow-up uses univariate ANOVAs with protection, discriminant function analysis to describe which combination of variables separates the groups, or simultaneous confidence intervals. Multivariate versions of every univariate design exist: randomised blocks, factorials and split plots, with a separate SSCP matrix for each source of variation.

MANOVA is preferred to separate ANOVAs when the variables are correlated and interest lies in the joint effect, for example the effect of an irrigation regime on a set of quality traits. Running p separate ANOVAs at level α inflates the family-wise error rate and ignores correlations that may make a combination of variables highly discriminating even when no variable is on its own.

Remember

6. Principal component analysis

Definition and computation

Principal component analysis (PCA) replaces p correlated variables by p uncorrelated linear combinations, the principal components, ordered by decreasing variance. The first component Y₁ = a₁′X maximises Var(a′X) = a′Σa subject to a′a = 1; the second maximises variance subject to being uncorrelated with the first, and so on. The solution: the coefficient vectors aⱼ are the eigenvectors of Σ, and Var(Yⱼ) = λⱼ, the corresponding eigenvalue, with λ₁ ≥ λ₂ ≥ ... ≥ λₚ. The components are uncorrelated, and the total variance is preserved: Σλⱼ = tr(Σ). The proportion of total variance explained by the j-th component is λⱼ/Σλₖ.

PCA can be based on the covariance matrix or the correlation matrix. Covariance-based PCA is dominated by variables with large variances, which depends on their units; correlation-based PCA (on standardised variables) treats all variables equally, and its eigenvalues sum to p. Agricultural data with mixed units (yield, height, days) are normally analysed on the correlation scale. For a 2 × 2 correlation matrix with correlation r, the eigenvalues are 1 + r and 1 − r, with eigenvectors (1, 1)/√2 and (1, −1)/√2: the first component is the sum and the second the contrast of the standardised variables.

The loadings, often reported as correlations between variables and components (aⱼₖ√λⱼ/σₖ for the covariance case, aⱼₖ√λⱼ for the correlation case), show which variables each component represents. Scores Yⱼ = aⱼ′(x − x̄) for each unit can be plotted, and biplots display units and variables together.

How many components, and uses

Common rules for retaining components are: enough to explain a stated proportion of variance (70-90%, a subjective choice); the Kaiser criterion for correlation-matrix PCA, keeping components with eigenvalues greater than 1 (the average); and the scree plot, retaining components before the 'elbow' where eigenvalues level off. Parallel analysis compares eigenvalues with those from random data of the same size and is more reliable. Example: with correlation-matrix eigenvalues 2.6, 1.2, 0.6, 0.4 and 0.2 from five correlated yield components (summing to 5), the first two explain (2.6 + 1.2)/5 = 76% of the variance, and the Kaiser rule also keeps two.

Uses in agriculture: summarising many morphological traits of germplasm accessions into two or three components for plotting and grouping; constructing indices from correlated soil or climate variables; reducing collinear predictors before regression (principal component regression); and analysing genotype-by-environment interaction (the AMMI model applies PCA to the interaction matrix, and GGE biplots to genotype plus interaction effects). PCA is a descriptive tool: it needs no distributional assumption, but its components need not have any causal or biological meaning, and the leading components capture variance, not necessarily the differences that matter for a particular purpose (that is the task of discriminant analysis).

Remember

7. Factor analysis

The factor model

Factor analysis explains the correlations among p observed variables by a smaller number m of unobserved common factors. The orthogonal factor model writes X − μ = LF + ε, with F the m common factors (mean 0, dispersion I), L the p × m matrix of factor loadings and ε the specific factors (mean 0, diagonal dispersion Ψ), independent of F. It implies the covariance structure Σ = LL′ + Ψ. For each variable, Var(Xⱼ) = hⱼ² + ψⱼ, where the communality hⱼ² = Σₖ lⱼₖ² is the part of the variance shared with other variables and ψⱼ the specific (unique) variance. For standardised variables, hⱼ² + ψⱼ = 1: a variable with loadings 0.8 and 0.3 on two factors has communality 0.64 + 0.09 = 0.73 and uniqueness 0.27.

The loadings are not unique: for any orthogonal matrix T, LT gives the same LL′. This rotational indeterminacy is used to seek interpretable solutions: varimax rotation (orthogonal) maximises the variance of squared loadings within factors, making each factor load heavily on a few variables; oblique rotations (promax, oblimin) allow correlated factors.

Estimation, choice of m and contrast with PCA

Loadings are estimated by the principal component (or principal factor) method, which takes the leading eigenvectors of S or R scaled by the square roots of their eigenvalues, or by maximum likelihood under normality. Maximum likelihood gives a large-sample test of the adequacy of m factors, with ((p − m)² − (p + m))/2 degrees of freedom; for p = 6 variables and m = 2 factors this is (16 − 8)/2 = 4. The number of factors is also chosen by eigenvalue rules and interpretability. Heywood cases, in which an estimated uniqueness is zero or negative, signal problems such as too many factors or too small a sample. Factor scores for individuals are estimated by the regression (Thomson) or weighted least squares (Bartlett) method.

Factor analysis differs from PCA in purpose and model. PCA is a transformation of the variables that explains total variance, with no error model; factor analysis is a statistical model that explains covariance (shared variance) through latent factors, setting aside specific variance. Components are exact linear functions of the observed variables; factors are not. In agriculture, factor analysis is used to identify underlying 'yield-component' or 'quality' factors among many correlated traits, to construct indices of soil fertility or of farm-household characteristics from survey items, and, in its confirmatory form (structural equation models), to test hypothesised relations among latent constructs such as adoption intention in extension research.

Remember

8. Canonical correlation analysis

Definition

Canonical correlation analysis studies the relationship between two sets of variables, X₁ (p variables) and X₂ (q variables), measured on the same units: for example, a set of soil properties and a set of crop responses, or yield components and quality traits. It finds linear combinations U₁ = a₁′X₁ and V₁ = b₁′X₂ with the largest possible correlation ρ₁, the first canonical correlation; then a second pair (U₂, V₂), uncorrelated with the first pair, with the next largest correlation, and so on, up to min(p, q) pairs.

The squared canonical correlations ρᵢ² are the eigenvalues of Σ₁₁⁻¹Σ₁₂Σ₂₂⁻¹Σ₂₁ (equivalently of Σ₂₂⁻¹Σ₂₁Σ₁₁⁻¹Σ₁₂), and the coefficient vectors are the corresponding eigenvectors. Canonical correlations are invariant to non-singular linear transformations within each set, in particular to changes of units. The first canonical correlation is at least as large as the absolute value of any simple correlation between a variable in one set and a variable in the other, and at least as large as any multiple correlation of one variable with the other set. Multiple correlation is the special case q = 1, and simple correlation p = q = 1.

Testing and interpretation

Bartlett's test of whether all canonical correlations are zero uses −[n − 1 − (p + q + 1)/2] Σ ln(1 − rᵢ²), approximately χ² with pq degrees of freedom; sequential versions test whether the remaining correlations after the first k are zero, with (p − k)(q − k) df. The test of all ρᵢ = 0 is equivalent to the test of independence of the two sets, which is a MANOVA-type (Wilks) test of the regression of one set on the other.

Interpretation uses structure coefficients (correlations between each variable and the canonical variates), which are more stable than the raw coefficients when variables within a set are correlated, and redundancy indices, which measure how much of one set's variance is explained by the other set's canonical variates. Large canonical correlations can arise from a single pair of variables and do not guarantee that the sets share much variance; redundancy analysis makes this explicit. Canonical correlation needs n well above p + q for stable results, and with many variables it overfits, which is why regularised versions are used for high-dimensional data such as linking genotype and phenotype sets.

Remember

9. Discriminant analysis and classification

Fisher's linear discriminant function

Discriminant analysis builds a rule for assigning a unit to one of several known groups from its measurements. Fisher's approach seeks the linear combination y = a′x that maximises the ratio of between-group to within-group variance. For two groups with a common dispersion matrix, the solution is a = S_p⁻¹(x̄₁ − x̄₂), the linear discriminant function (LDF). A new unit x is assigned to group 1 if a′x exceeds the midpoint m = (1/2)a′(x̄₁ + x̄₂), and to group 2 otherwise, when prior probabilities and misclassification costs are equal.

The separation achieved is measured by the Mahalanobis distance: the variance of y is a′S_pa = D², and the difference between the group means of y is also D² = (x̄₁ − x̄₂)′S_p⁻¹(x̄₁ − x̄₂). For two multivariate normal populations with common Σ and population distance Δ, the optimal rule misclassifies with probability Φ(−Δ/2) in each direction (equal priors); with Δ = 2 this is Φ(−1) = 0.159, and with Δ = 4 it is Φ(−2) = 0.023. Seed-lot identification, separating varieties by morphological or spectral measurements, classifying soils into fertility classes and identifying diseased plants from reflectance are typical applications.

Bayes rules, unequal covariances and error rates

The general rule minimises expected cost of misclassification: assign x to group 1 if f₁(x)/f₂(x) > [c(1 | 2)p₂]/[c(2 | 1)p₁], where pᵢ are prior probabilities and c(i | j) the cost of assigning a group-j unit to group i. With normal populations and a common Σ this gives a linear rule (LDA) with a shifted cut-off; with unequal covariance matrices it gives quadratic discriminant analysis (QDA), with quadratic boundaries and more parameters to estimate. For g groups the rule assigns to the group with the largest posterior probability, and Fisher's canonical discriminant functions (eigenvectors of W⁻¹B, at most min(p, g − 1) of them) give low-dimensional plots of group separation.

The apparent error rate, obtained by classifying the same units used to build the rule, is optimistic (biased downward). Honest estimates use a separate test sample or cross-validation, of which Lachenbruch's holdout (leave-one-out) method is the classic version. Logistic regression is an alternative to LDA that needs no normality of the predictors and handles categorical predictors; LDA is more efficient when normality holds. Modern classifiers (classification trees, random forests, support vector machines) are used for high-dimensional data such as spectral or image features, and are judged by the same cross-validated error rates.

Remember

10. Cluster analysis

Distances and hierarchical clustering

Cluster analysis groups units that are similar on many variables when no groups are known in advance. It starts from a measure of dissimilarity: Euclidean distance (after standardising variables measured in different units), Mahalanobis distance, Manhattan distance, or, for binary or categorical data such as molecular-marker profiles, similarity coefficients (simple matching, Jaccard) converted to distances.

Agglomerative hierarchical clustering starts with every unit in its own cluster and repeatedly merges the two closest clusters. Methods differ in how the distance between clusters is defined: single linkage (nearest neighbour) uses the smallest distance between members and tends to produce long chains; complete linkage (farthest neighbour) uses the largest distance and produces compact clusters; average linkage (UPGMA) uses the mean distance between all pairs and is widely used for marker-based diversity trees; Ward's method merges the pair that gives the smallest increase in the total within-cluster sum of squares and tends to produce clusters of similar size. The result is displayed as a dendrogram, cut at a chosen height to form groups. The cophenetic correlation, between the original distances and the distances implied by the dendrogram, measures how faithfully the tree represents the data.

k-means and practical choices

k-means clustering partitions the units into a pre-specified number k of clusters so as to minimise the total within-cluster sum of squared Euclidean distances to the cluster means. The algorithm alternates between assigning each unit to the nearest centre and recomputing the centres; it converges to a local optimum that depends on the starting centres, so several random starts are used. Choosing k uses the elbow in the within-cluster sum of squares, the average silhouette width or the gap statistic. Model-based clustering (finite mixtures of multivariate normals) provides likelihood-based choice of k (for example by BIC) and probabilities of membership.

Results depend on the variables chosen, their scaling, the distance and the method, so cluster analysis is exploratory: different reasonable choices should give broadly similar groupings before the groups are interpreted. A small example of single versus complete linkage: with distances d(A, B) = 2, d(A, C) = 6, d(B, C) = 5, after merging A and B, the distance from cluster AB to C is 5 under single linkage, 6 under complete linkage and 5.5 under average linkage. Clustering of germplasm accessions, of farms into typologies for targeting extension, and of environments into mega-environments for variety testing are common agricultural applications.

Remember

11. Mahalanobis D² and genetic divergence

D² as a generalised distance

The Mahalanobis distance between two populations with means μ₁, μ₂ and common dispersion Σ is Δ² = (μ₁ − μ₂)′Σ⁻¹(μ₁ − μ₂), estimated by D² = (x̄₁ − x̄₂)′S⁻¹(x̄₁ − x̄₂) with S the pooled within-population dispersion matrix. Unlike Euclidean distance it accounts for the different variances of the variables and for their correlations: two correlated traits that both differ between populations are not counted twice. D² is invariant to non-singular linear transformations of the variables, so units of measurement do not matter. For uncorrelated variables it reduces to the sum of squared standardised differences, Σ(δⱼ/σⱼ)². With two traits of unit variance and correlation r, D² = (δ₁² − 2rδ₁δ₂ + δ₂²)/(1 − r²); for differences of 1 and 1 with r = 0.5 this gives (1 − 1 + 1)/0.75 = 1.33, compared with 2 if the correlation were ignored.

Its significance for two groups is tested through the two-sample Hotelling T² = [n₁n₂/(n₁ + n₂)]D² (Chapter 4); with many populations, Rao's procedures test the significance of each pairwise D² and of the overall divergence, usually through Wilks' Λ.

Genetic divergence analysis in plant breeding

C. R. Rao's use of D² for classifying populations became a standard tool in Indian and South Asian plant breeding for assessing genetic diversity among genotypes before choosing parents for hybridisation. The usual procedure is: measure several quantitative traits on all genotypes in a replicated trial; test for differences by Wilks' Λ; transform the correlated traits to uncorrelated ones (pivotal condensation of the pooled error dispersion matrix) so that D² is a sum of squares; compute D² for all pairs; group the genotypes into clusters, classically by Tocher's method, which adds genotypes to a cluster so long as the average intra-cluster distance stays below the largest of the minimum distances; and compute average intra- and inter-cluster distances.

Breeders then choose parents from clusters separated by large inter-cluster distances, expecting more heterosis and wider segregation, while also considering cluster means for the traits of interest. The relative contribution of each trait to divergence is often summarised by how frequently it ranks first among the pairwise differences. Two cautions apply: genetic divergence measured on phenotypes in one environment includes environmental and genotype-by-environment effects, and geographical diversity does not necessarily imply genetic diversity, a point repeatedly found in such studies. Molecular-marker distances (Nei's, Jaccard-based) now complement morphological D², and the two often give different groupings, since they reflect different parts of the genome and different evolutionary histories.

Remember

12. Multidimensional scaling and correspondence analysis

Multidimensional scaling

Multidimensional scaling (MDS) represents n units as points in a low-dimensional space so that the distances between points approximate a given matrix of dissimilarities. Classical (metric, Torgerson) MDS starts from the matrix of squared distances D², forms B = −(1/2)JD²J, where J = I − (1/n)11′ is the centring matrix, and uses the leading eigenvectors of B scaled by the square roots of their eigenvalues as coordinates. When the dissimilarities are Euclidean distances computed from a data matrix, classical MDS gives the same configuration as PCA of that data (principal coordinate analysis is another name for classical MDS). Negative eigenvalues of B indicate that the dissimilarities are not exactly Euclidean.

Non-metric MDS (Kruskal) uses only the rank order of the dissimilarities and finds a configuration minimising stress, a measure of mismatch between the ordered dissimilarities and the fitted distances; stress below about 0.1 is commonly regarded as a fair fit and below 0.05 as good, though such guidelines are conventions. MDS is used to display diversity among accessions from marker-based distances, similarity among environments in variety trials, and perceptual similarity of products in sensory studies.

Correspondence analysis

Correspondence analysis (CA) is the analogue of PCA for a two-way contingency table of counts, such as crop varieties by agro-climatic zones in a germplasm survey, or farm type by adopted practices. It represents the row and column categories as points in a low-dimensional space in which proximity reflects similarity of profiles. Row profiles (each row divided by its total) are compared with chi-square distances, which weight each column by the inverse of its overall proportion.

The total inertia of the table equals the Pearson chi-square statistic divided by the grand total, χ²/n, and the analysis decomposes it into principal inertias (squared singular values of the matrix of standardised residuals), the CA counterpart of eigenvalues. The proportion of inertia shown by the first two dimensions measures the quality of the map. Example: a table with n = 500 and χ² = 60 has total inertia 0.12; if the first two principal inertias are 0.07 and 0.03, the map displays 83% of the inertia. Rows and columns can be plotted together (a symmetric map), but distances between a row point and a column point are not directly interpretable; association is read from directions and from the positions relative to the origin. Multiple correspondence analysis extends CA to several categorical variables, as in farm-household surveys.

Remember

What the exams ask

How it is used, across a career

As a student

A student must handle the matrix algebra and the main results: conditional distributions of the MVN, the F transformation of T², Wilks' Λ and its approximation, eigen-analysis for PCA and canonical correlations, communalities in factor analysis, the Fisher discriminant and its error rate, and D² with Tocher clustering. AICE-JRF names most of these topics; plant-breeding theses routinely use D², PCA and cluster analysis, and need them done and interpreted correctly.

Five years in

A statistician five years into work at a breeding or soil-science station runs multivariate analyses every season: MANOVA on sets of quality traits, PCA and biplots of germplasm and of G × E tables, D² divergence and clustering for parent selection, discriminant rules for identifying varieties or seed lots from spectral data with cross-validated error rates, and factor or correspondence analysis of farm-survey data. Checking scaling, outliers, sample size relative to the number of variables, and the stability of clusters is part of the job.

Twenty-five years in

A senior statistician with twenty-five years of experience sets standards for multivariate work in a programme: which variables and scalings are used in diversity studies, how clusters and components are validated before being reported, how classification rules for seed certification or disease detection are tested before adoption, and how high-dimensional data (genomic, spectral, remote sensing) are handled with regularised methods. Such a statistician reviews multivariate claims in breeding and survey reports and guards against over-interpretation of components and clusters.

PhD, twenty-five years in

A PhD with long experience develops multivariate methods for agricultural data: high-dimensional covariance estimation for genomic prediction, multivariate mixed models for correlated traits across environments, model-based clustering and mixture models for germplasm, sparse PCA and canonical correlation for linking omics layers, and the theory of test statistics when p is comparable to n. Such a researcher proves properties, checks them by simulation and supervises doctoral students.

Written from

Derivations of the Wishart density and of the exact distributions of MANOVA statistics are not given. Not covered in depth: structural equation models, partial least squares, independent component analysis, high-dimensional (p > n) theory, robust multivariate methods, and machine-learning classifiers beyond a mention. Rules of thumb for retaining components, for stress values in MDS and for cluster numbers are conventions that vary between authors. The relative-contribution method for traits in D² studies is described in general terms because published variants differ. Statistical tables are not reproduced; questions state the values they need.

Original practice questions written to the syllabus and pattern of the exams named in each question's examTags. They are not past papers and do not reproduce any board's questions. Every numerical answer was computed with a script; the values needed are stated in the question.

Practise it, and prove it

These notes are free to read here, with no account.

What the Laisha app adds:

Free for ten days, then a subscription. We would rather write that here than have you find it later.

Open Laisha →