Why Study Multivariate Data?

• In ecology, often measure many variables on the same organisms:
  • Plant morphology: height, width, leaf size, flower size
  • Environmental conditions: temperature, pH, nutrients
  • Species composition: presence/absence, abundance
• Challenge: How do we make sense of all these variables at once?
• PCA helps by:
  • Reducing many correlated variables to a few key dimensions
  • Visualizing patterns in high-dimensional data
  • Identifying the most important sources of variation

Darlingtonia californica (Cobra Lily)

The Cobra Lily Dataset

• Analyze morphological measurements from 89 individual Darlingtonia californica plants collected from 4 sites in the Siskiyou Mountains of Oregon.
• Variables measured (10 total):
  • Morphological measurements (mm):
    • height - Height of pitcher
    • mouth_diam - Mouth diameter
    • tube_diam - Tube diameter
    • keel_diam - Keel diameter
    • wing1_length, wing2_length
    • Wing lengths (2 measurements)
    • wingsprea - Wing spread
  • Biomass measurements (g):
    • hoodmass_g - Hood mass
    • tubemass_g - Tube mass
    • wingmass_g - Wing mass

![Darlingtonia pitcher structure showing measured features\]

![](images/clipboard-4005369252.png){width="207"}

- **Sites**:
-   TJH: T.J. Howell's Fen
-   DG: Day's Gulch
-   LEH: L.E. Hunter's Fen
-   HD: High Divide
 

Why Standardization Matters

• Without standardization:
  • Variables with larger values (like height in mm) would dominate the analysis over variables with smaller values (like biomass in g).
• Example from data:
  • Height: mean ≈ 600 mm, range 322-845 mm -
  • Wing mass: mean ≈ 0.16 g, range 0.02-5.00 g
• calculated distances without standardizing, height would have 1000× more influence than wing mass!
• Solution: Z-score standardization
  • \(Z = \frac{(Y_i - \bar{Y})}{s}\)
• transforms each variable to:
  • Mean = 0
  • Standard deviation = 1
  • Units: “standard deviations from the mean”

• Before standardization: - Height dominates (larger scale)
• After standardization:
  • All variables contribute equally
  • Measured in “standard deviations”

# Example: standardize height
darlingtonia %>%
  summarize(
    mean_height = mean(height),
    sd_height = sd(height),
    mean_wingmass = mean(wingmass_g),
    sd_wingmass = sd(wingmass_g)
  )
## # A tibble: 1 × 4
##   mean_height sd_height mean_wingmass sd_wingmass
##         <dbl>     <dbl>         <dbl>       <dbl>
## 1        616.      101.         0.154       0.102

# After standardization, both will have:
# mean = 0, sd = 1

Euclidean Distance Between Individuals

• How different are two plants?
  • For two plants measured on multiple variables, calculate Euclidean distance:

    • For 2 variables (height and spread):

      • \(d_{i,j} = \sqrt{(y_{i,1} - y_{j,1})^2 + (y_{i,2} - y_{j,2})^2}\)

    • For n variables:

      • \(d_{i,j} = \sqrt{\sum_{k=1}^{n}(y_{i,k} - y_{j,k})^2}\)

    • Pythagorean theorem extended to many dimensions!

    • Key point: Distance should use standardized values so all variables contribute equally

Step 2: Z-score Transformation

Visualizing Standardized Data

Why Correlations Matter for PCA

• PCA works best when variables are correlated:
  • If variables are uncorrelated → PCA provides no benefit
  • If variables are highly correlated → PCA can reduce dimensions effectively
• In data, expect correlations because:
  • Larger plants have larger measurements overall
  • Wing dimensions likely correlate with each other
  • Mass measurements correlate with size measurements
• The correlation matrix shows relationships between all pairs of variables.
• Strong correlations (positive or negative) suggest PCA will be useful!

Strong positive correlations (dark blue) suggest variables measure similar aspects of plant size/shape.

The Concept of PCA

• PCA creates new variables (components) that are:
  • Linear combinations of original variables
  • Uncorrelated with each other
  • Ordered by importance (amount of variance explained)
• The mathematical idea:
  • For each plant i and component k:
  • \(Z_{ik} = a_{i1}Y_1 + a_{i2}Y_2 + ... + a_{in}Y_n\)
• Where:
  • \(Z_{ik}\) = score on principal component k
  • \(Y_j\) = standardized value of original variable j
  • \(a_{ij}\) = loading (weight) of variable j on component i
• Goal: Find loadings that maximize variance in first component, then second, etc.

Conceptual diagram of PCA: reducing 2 variables to 1 principal axis

• The first principal axis passes through maximum variance.
• Key properties: - As many components as variables - First components explain most variation - Later components explain residual variation

Understanding Components are:

• 1. Orthogonal (perpendicular)
  • Components are uncorrelated
  • Each captures independent information
• 2. Ordered by variance
  • PC1 explains the most variance
  • PC2 explains second most (of remaining variance)
  • And so on…
• 3. Linear combinations
  • Each component is a weighted sum of original variables
  • Weights (loadings) show variable importance

Why this matters:

• Dimension reduction:
  • Start with 10 correlated variables
  • End with 2-3 uncorrelated components
  • Retain most information with fewer dimensions
• Easier interpretation:
  • PC1 might represent “overall size”
  • PC2 might represent “shape”
  • Components often have biological meaning
• Better for analysis:
  • Can use in regression, ANOVA
  • No multicollinearity problems
  • Variables are standardized

How Much Variance Does Each PC Explain?

• Eigenvalues (\(\lambda\)) measure variance explained by each component.
• The proportion of variance for PC j is:
  • \(\text{Proportion}_j = \frac{\lambda_j}{\sum_{k=1}^{n}\lambda_k}\)
• In practice:

  • PC1 always explains the most variance

  • We want the first few PCs to explain most (>70%) of total variance

  • allows dimension reduction

How Many Components to Use?

• A scree plot shows variance explained by each component.
• How to read it:
  • Look for an “elbow” - sharp change in slope
  • Keep components before the elbow
  • Discard components in the “scree” (rubble)
• In this data:
  • PC1 explains 49.9% of variance
  • PC2 explains 17.7%
  • PC3 explains 15.3%
  • First 3 PCs explain 77.3% of variance
• Decision: Use PC1, PC2, and PC3 for interpretation.
• remaining 7 components explain 17.1% (mostly noise).

What are Loadings?

• Loadings are the weights (\(a_{ij}\)) that show how each original variable contributes to each principal component.
• For PC1:
  • \(Z_1 = a_{11}Y_1 + a_{12}Y_2 + ... + a_{1p}Y_p\)
• Interpreting loadings:
  • Large positive loading: variable increases with PC
  • Large negative loading: variable decreases with PC
  • Near-zero loading: variable unimportant for PC
• Guidelines:
  • Loadings > 0.3 or < -0.3 are worth considerin
  • Loadings > 0.5 or < -0.5 are important
  • Sign (+ or -) can flip; patterns matter more

• PC1 loadings:
  • Most variables have similar positive values
  • All ~0.31 to 0.40 except keel (-0.18) and tube (-0.002)
  • Suggests PC1 measures “overall size”

##                 PC1    PC2    PC3
## height       -0.278 -0.449  0.198
## mouth_diam   -0.369 -0.110  0.214
## tube_diam     0.040 -0.602 -0.390
## keel_diam     0.202 -0.491 -0.328
## wing1_length -0.375  0.171 -0.280
## wing2_length -0.342  0.131 -0.415
## wingsprea    -0.239  0.156 -0.549
## hoodmass_g   -0.385 -0.057  0.185
## tubemass_g   -0.355 -0.320  0.263
## wingmass_g   -0.393  0.080 -0.007

PC2 Interpretation

• Looking at PC2 loadings:
• Key patterns:
  • Height: -0.419 (negative, moderate)
  • Mouth: -0.250 (negative, moderate)
  • Wing measurements: positive (0.275-0.434)
  • Wing spread: 0.434 (positive, strong)
  • Biomass: mostly negative
• Biological meaning: PC2 represents shape trade-offs
  • tall narrow pitchers vs. short wide pitchers with large wings.

## # A tibble: 10 × 2
##    Variable     Loading
##    <chr>          <dbl>
##  1 tube_diam    -0.602 
##  2 keel_diam    -0.491 
##  3 height       -0.449 
##  4 tubemass_g   -0.320 
##  5 wing1_length  0.171 
##  6 wingsprea     0.156 
##  7 wing2_length  0.131 
##  8 mouth_diam   -0.110 
##  9 wingmass_g    0.0800
## 10 hoodmass_g   -0.0575

• High PC2 → short plants with large wings
• Low PC2 → tall narrow plants

PC3 Interpretation

• Looking at PC3 loadings:
• Key patterns:
• Dominated by:
  • Tube diameter: 0.743 (very strong positive)
  • Keel diameter: 0.576 (strong positive)
  • Most other variables near zero
• Biological meaning: PC3 represents tube girth
  • width/diameter of tube structure independent of height or wing size.
  • Plants with high PC3 scores = fat tubes, allowing to trap larger prey?

## # A tibble: 10 × 2
##    Variable      Loading
##    <chr>           <dbl>
##  1 wingsprea    -0.549  
##  2 wing2_length -0.415  
##  3 tube_diam    -0.390  
##  4 keel_diam    -0.328  
##  5 wing1_length -0.280  
##  6 tubemass_g    0.263  
##  7 mouth_diam    0.214  
##  8 height        0.198  
##  9 hoodmass_g    0.185  
## 10 wingmass_g   -0.00743

• High PC3 → fat tubes
• Low PC3 → skinny tubes

What are PC Scores?

• PC scores are values of each principal component for each obs./plant
• Each plant gets a score on each PC:
  • PC1 score = measure of size
  • PC2 score = measure of shape
  • PC3 score = measure of tube girth
• How they’re calculated:
• For plant i on PC1:
  • \(\text{PC1}_i = 0.31×height_i + 0.40×mouth_i + ...\)
  • Using the standardized values and loadings we saw earlier.

Ordination of Sites

• We can plot PC scores to visualize:
  • How plants are distributed in PC space
  • Whether sites differ in morphology
  • Patterns and groupings
• In this plot:
  • Each point = one plant
  • Position determined by PC1 and PC2 scores
  • Colors indicate collection site
• Observations:
  • Sites show some separation
  • HD (orange) plants have lower PC1 (smaller)
  • TJH (pink) plants variable but lower PC2
  • Some overlap among sites
• This is called ordination - ordering objects along axes.

The ordination

Reading a Biplot

• A biplot shows both:
  • Points = plants (PC scores)
  • Arrows = original variables (loadings)
• How to interpret arrows:
• Length = how important the variable is
• Direction = how the variable relates to PCs
• Angle between arrows = correlation
  • Small angle = positively correlated
  • 90° = uncorrelated
  • 180° = negatively correlated
• How to interpret points:
  • Plants in the direction of an arrow are high in that variable
  • Plants opposite an arrow are low in that variable

The plot

Key Patterns in the Biplot

• Variable relationships:
  • Most morphological and biomass variables point together (right and slightly down):

    • Height, mouth, masses - These are highly correlated

    • All contribute to PC1 (size)

• Wing variables (wing1_length, wing2_length, wingspread) point more upward:
  • Contribute positively to PC2
    • Represent the “wing dimension” of shape
• Site separation:

  • HD (blue): Lower PC1, plants are smaller overall

  • DG (orange): Higher PC1, plants are larger

  • TJH (yellow) and LEH (green): Intermediate and variable

Biological Interpretation

• What determines plant morphology?
  • Size (PC1, 45.8% variance)

    • Primary source of variation

    • Most measurements scale together

    • May reflect resource availability or age

  • Shape trade-offs (PC2, 16.7% variance)
    • Tall narrow vs. short wide with large wings
    • May reflect different trapping strategies
    • Or developmental/environmental constraints
  • Tube dimensions (PC3, 14.8% variance)
    • Independent aspect of morphology
    • May affect prey capture efficiency
  • Sites differ primarily in overall size, with some shape variation.

Box Plots of PC Scores

HD plants are significantly smaller than plants at other sites.

• TJH plants have different shapes

  • taller and narrower with smaller wings.

Calculating Euclidean Distance Between Groups

• Calculate the distance between site centroids (mean PC scores):
  • \(d_{k,j} = \sqrt{\sum_{i=1}^{p}(\bar{Z}_{i,k} - \bar{Z}_{i,j})^2}\)
  • Where \(\bar{Z}_{i,k}\) is the mean of PC i for site k.
• Using all PCs preserves the true Euclidean distance from the original data.
• Using just PC1-PC3 (77% of variance) gives approximate distances.

HD is most different from other sites (smallest plants).

# Calculate site centroids (means)
site_means <- pc_scores %>%
  group_by(site) %>%
  summarize(across(PC1:PC3, mean))

site_means
## # A tibble: 4 × 4
##   site     PC1    PC2    PC3
##   <chr>  <dbl>  <dbl>  <dbl>
## 1 DG    -1.02   0.604  0.542
## 2 HD     2.63  -0.863 -1.26 
## 3 LEH   -0.659  0.155 -0.555
## 4 TJH    0.423 -0.345  0.618

# Calculate distance matrix
site_dist <- site_means %>%
  dplyr::select(PC1, PC2, PC3) %>%
  dist() %>%
  round(3)
site_dist
##       1     2     3
## 2 4.328            
## 3 1.241 3.511      
## 4 1.732 2.942 1.672

What We Learned

• Dimension Reduction Success
  • Started with: 10 correlated variables
    
  • Reduced to: 3 principal components
    
  • Retained: 77.3% of variance
• This makes data easier to visualize, interpret, and analyze.
• Biological Insights - Three major axes of morphological variation: - PC1 (45.8%): Overall size
  • Larger vs. smaller plants
  • HD site has smallest plants - PC2 (16.7%): Shape contrast
  • Tall/narrow vs. short/wide with large wings
  • TJH site has taller, narrower plants - PC3 (14.8%): Tube girth
  • Fat vs. skinny tubes
  • No significant site differences

• Statistical Findings
  • Sites differ significantly in size (PC1)
  • Sites differ significantly in shape (PC2)
  • Sites do not differ in tube girth (PC3)
• PCA Advantages
  • Reduced 10 dimensions to 3
    
  • Created uncorrelated variables
    
  • Revealed biological patterns
    
  • Enabled statistical testing
    
  • Simplified visualization
• Key Take-Home Messages
  • PCA good for exploratory analysis of multivariate data
  • Always standardize when different units/scales
  • Loadings reveal variables contributions to each component
  • Scores used in subsequent analyses (ANOVA, regression)
  • Biplots visualize both samples and variables simultaneously

• Functional PCA: - For time series or spatial data - Treats curves as observations - Common in physiology
• Probabilistic PCA: - Includes uncertainty estimation - Can handle missing data - Maximum likelihood framework
• Kernel PCA: - Handles non-linear relationships - Uses kernel trick from machine learning - More flexible but harder to interpret
  
• These are beyond our scope but worth knowing exist!

❌ Using PCA on categorical data
→ ✓ Use appropriate method (CA, MCA)

❌ Not checking for outliers
→ ✓ Plot scores, check for extreme values

❌ Assuming PCs have biological meaning
→ ✓ Interpret carefully, components are mathematical

❌ Using PCA with too few samples
→ ✓ Need n > 3p (preferably n > 50)

❌ Not reporting variance explained
→ ✓ Always report % variance for each PC

❌ Forgetting that PCA is exploratory
→ ✓ Use for patterns, not hypothesis testing