High-Dimensional Data

Lecture 7

Shannon Gallagher

September 16, 2026

Code 
source("_setup.R")
metrics <- state_metrics |>
  transmute(
    state = state_upper,
    event_records,
    injuries = log1p(injuries),
    deaths = log1p(deaths),
    damage = log1p(property_damage),
    duration = mean_duration
  )
x <- metrics |> select(-state) |> scale()

Reminders, previously, and today…

  • explain why distance depends on scale and variables,
  • read a multidimensional-scaling display,
  • understand PCA as a rotated coordinate system,
  • use loadings and explained variance to interpret components.

Thinking about distance…

Euclidean distance between observations \(i\) and \(j\) is:

\[ d(i,j)=\sqrt{\sum_{k=1}^{p}(x_{ik}-x_{jk})^2} \]

The formula is easy. Choosing variables and scales is the hard part.

Distances in general

  • Euclidean: straight-line separation; sensitive to scale.
  • Manhattan: total coordinate-wise separation; less dominated by one large step.
  • Correlation distance: similarity of profiles after ignoring overall level.
  • Domain-specific distance: weights or transforms chosen for the actual decision.

There is no universal distance for “similar states.”

Inspect the distance matrix

Code 
distance_matrix <- as.matrix(dist(x))
round(distance_matrix[1:6, 1:6], 2)
     1    2    3    4    5    6
1 0.00 1.80 3.05 1.03 2.73 1.53
2 1.80 0.00 3.49 1.96 4.19 1.77
3 3.05 3.49 0.00 2.55 2.76 3.05
4 1.03 1.96 2.55 0.00 2.35 1.67
5 2.73 4.19 2.76 2.35 0.00 3.10
6 1.53 1.77 3.05 1.67 3.10 0.00

The matrix is symmetric, has zero on the diagonal, and contains one distance for every pair of observations.

Dollars dominate without scaling

If we combine raw event counts, deaths, hours, and billions of dollars, damage controls almost every distance.

Standardization makes one standard deviation count equally in each variable—but that is also a substantive choice.

Checkpoint 1 · Define “similar”

Two states have similar total damage but very different deaths and event counts. Are they similar? Choose the variables and scaling that match one plausible client question.

Multi-dimensional scaling (MDS)

MDS starts from pairwise distances and finds a low-dimensional arrangement that preserves those distances as closely as possible.

Nearby points are similar under the chosen distance. Axis directions usually have no standalone meaning.

MDS example with NOAA state profiles

Code 
mds <- cmdscale(dist(x), k = 2)
mds_df <- tibble(state = metrics$state, dim1 = mds[, 1], dim2 = mds[, 2])

ggplot(mds_df, aes(dim1, dim2, label = state)) +
  geom_hline(yintercept = 0, color = "gray85") +
  geom_vline(xintercept = 0, color = "gray85") +
  geom_text(color = deep_gold, check_overlap = TRUE) +
  labs(title = "Nearby states have similar standardized NOAA profiles", x = "MDS dimension 1", y = "MDS dimension 2")
Multidimensional scaling plot of state-level NOAA profiles.

View structure with additional variables

Code 
raw_x <- metrics |> select(-state)
raw_mds <- cmdscale(dist(raw_x), k = 2)

ggplot(tibble(state = metrics$state, x = raw_mds[, 1], y = raw_mds[, 2]), aes(x, y, label = state)) +
  geom_text(color = blue_gray, check_overlap = TRUE) +
  labs(title = "Raw units", x = NULL, y = NULL)
ggplot(mds_df, aes(dim1, dim2, label = state)) +
  geom_text(color = deep_gold, check_overlap = TRUE) +
  labs(title = "Standardized variables", x = NULL, y = NULL)

Two MDS plots comparing raw and standardized state-level NOAA metrics.

Two MDS plots comparing raw and standardized state-level NOAA metrics.

Evaluate how well two dimensions preserve distance

Code 
embedded_distance <- as.vector(dist(mds))
original_distance <- as.vector(dist(x))

tibble(original_distance, embedded_distance) |>
  ggplot(aes(original_distance, embedded_distance)) +
  geom_point(alpha = 0.5, color = deep_gold) +
  geom_abline(slope = 1, intercept = 0, color = blue_gray) +
  coord_equal() +
  labs(title = "Departures from the diagonal show distortion", x = "Original standardized distance", y = "Distance in the 2D map")
Scatterplot of original distances and distances in the two-dimensional MDS solution.

What MDS does not say

  • It does not identify causal groups.
  • The axes may rotate or flip without changing the solution.
  • A two-dimensional view necessarily loses information.
  • Clusters can depend strongly on distance and preprocessing.

Dimension reduction—searching for variance

Principal component analysis finds directions that:

  1. are linear combinations of the original variables,
  2. capture as much variance as possible,
  3. are mutually orthogonal.

PCA rotates the coordinate system; it does not discover truth.

Principal Component Analysis (PCA)

For centered variables, the first component has the form

\[ PC_1=a_{11}X_1+a_{21}X_2+\cdots+a_{p1}X_p. \]

The weights are chosen so that the scores on \(PC_1\) have maximum variance, subject to the loading vector having length one. Later components maximize remaining variance while remaining orthogonal to earlier components.

What are principal components?

  • Centering moves every variable to mean zero.
  • Scaling gives every variable unit standard deviation.
  • PCA on a covariance matrix preserves original-unit variance.
  • PCA on a correlation matrix treats standardized variables equally.

With counts, deaths, hours, and dollars together, scaling is essential unless the units are intentionally weighted.

Computing Principal Components

Code 
pca <- prcomp(x, center = FALSE, scale. = FALSE)
summary(pca)
Importance of components:
                          PC1    PC2    PC3     PC4     PC5
Standard deviation     1.5895 1.0039 0.8378 0.70626 0.51494
Proportion of Variance 0.5053 0.2016 0.1404 0.09976 0.05303
Cumulative Proportion  0.5053 0.7068 0.8472 0.94697 1.00000

Making PCs interpretable with loadings

Code 
round(pca$rotation[, 1:3], 2)
               PC1   PC2   PC3
event_records 0.50  0.17 -0.58
injuries      0.55 -0.03 -0.29
deaths        0.43  0.23  0.69
damage        0.49 -0.07  0.30
duration      0.15 -0.95  0.05

A loading describes how strongly an original variable contributes to a component. Component signs can flip; relative patterns matter.

Making PCs interpretable with biplots

Code 
loading_scale <- 4
loading_df <- as_tibble(pca$rotation[, 1:2], rownames = "variable") |>
  mutate(PC1 = PC1 * loading_scale, PC2 = PC2 * loading_scale)
scores <- as_tibble(pca$x[, 1:2]) |>
  mutate(state = metrics$state)

ggplot(scores, aes(PC1, PC2)) +
  geom_text(aes(label = state), color = deep_gold, check_overlap = TRUE) +
  geom_segment(data = loading_df, aes(x = 0, y = 0, xend = PC1, yend = PC2),
               inherit.aes = FALSE, arrow = arrow(length = unit(0.15, "in")), color = blue_gray) +
  geom_text(data = loading_df, aes(PC1, PC2, label = variable), inherit.aes = FALSE,
            color = blue_gray, vjust = -0.5) +
  labs(title = "Scores locate states; arrows show loading directions")
PCA biplot of state scores and variable loading directions.

Arrows pointing together indicate variables that contribute similarly to these two components. Their displayed length depends on the biplot scaling convention.

Checkpoint 2 · Name PC1

Read the PC1 loadings. Propose a plain-language name for the component, then list one reason your name is incomplete or potentially misleading.

NOAA state profiles: PC1 and PC2

Code 
scores <- as_tibble(pca$x[, 1:2]) |>
  mutate(state = metrics$state)

ggplot(scores, aes(PC1, PC2, label = state)) +
  geom_hline(yintercept = 0, color = "gray85") +
  geom_vline(xintercept = 0, color = "gray85") +
  geom_text(color = deep_gold, check_overlap = TRUE) +
  labs(title = "PCA summarizes state profiles on new axes")
Scatterplot of state PCA scores on the first two components.

How many principal components should we use?

Code 
variance <- pca$sdev^2 / sum(pca$sdev^2)
tibble(component = seq_along(variance), variance) |>
  ggplot(aes(component, variance)) +
  geom_col(fill = gold) +
  geom_line() +
  geom_point() +
  scale_y_continuous(labels = percent) +
  labs(title = "A scree plot makes information loss visible", x = "Principal component", y = "Variance explained")
Scree plot showing variance explained by each principal component.

Create a scree plot (aka “elbow plot”) to choose

Code 
tibble(component = seq_along(variance), cumulative = cumsum(variance)) |>
  ggplot(aes(component, cumulative)) +
  geom_line(color = deep_gold, linewidth = 1) +
  geom_point(color = deep_gold, size = 3) +
  geom_hline(yintercept = 0.8, linetype = 2, color = blue_gray) +
  scale_y_continuous(labels = percent, limits = c(0, 1)) +
  labs(title = "Retaining components trades simplicity for information", x = "Number of components", y = "Cumulative variance explained")
Cumulative variance explained by the principal components.

PCA: singular value decomposition (SVD)

For a centered data matrix \(X\),

\[ X=UDV^T. \]

  • columns of \(V\) are the loading directions;
  • \(UD\) contains the component scores;
  • squared singular values determine variance explained.

This connection makes PCA computationally stable and links it to many other matrix methods.

Checkpoint 3 · Keep how many?

Using the scree plot, choose how many components you would retain for exploration. Report the cumulative variance and one important reason variance explained is not the only criterion.

Recap and next steps

  • Distance encodes a definition of similarity.
  • MDS preserves pairwise distances; its axes are usually not interpreted directly.
  • PCA preserves variance along linear directions.
  • Loadings, scores, and information loss must be read together.

Next: nonlinear embeddings and the transition to time.