Functional principal components analysis (FPCA)
Background on FPCA
FPCA characterizes modes of variability by decomposing functional observations into population level basis functions and subject-specific scores. The basis functions have a clear interpretation, analogous to that of PCA: the first basis function explains the largest direction of variation, and each subsequent basis function describes less. The FPCA model is typically written
\[ X_i(r) = \mu(r) + \sum_{k=1}^{K} c_{ik}\psi_{k}(r) + \epsilon_i(r) \]
where \(\mu(r)\) is the population mean, \(\psi_{k}(r)\) are a set of orthonormal population-level functional principal components, \(c_{ik}\) are subject-specific scores with mean zero and variance \(\lambda_k\), and \(\epsilon_i(r)\) are residual curves. Estimated functional principal components \(\widehat{\psi}_1(r), \widehat{\psi}_2(r), \ldots, \widehat{\psi}_{K}(r)\) and corresponding variances \(\widehat{\lambda}_1 \geq \widehat{\lambda}_2 \geq \ldots \geq \widehat{\lambda}_K\) are obtained from a truncated Karhunen-Lo`eve decomposition of the sample covariance \(\widehat{\Sigma}(s,t) = \widehat{\mbox{Cov}}(Y_i(s), X_i(r))\). The truncation lag \(K\) is often chosen so that the resulting approximation accounts for at least 95% of observed variance (Wrobel et al. (2016)).
Conceptually, functional principal components can be thought of as patterns in the data that explain the most variance, and scores \(c_{ik}\) are weights that describe how much each functional principal component contributes to the shape of a given subject’s spatial summary function.
Implementing FPCA
Load and visualize data
Here we load data from the VectraPolarisData ovarian cancer dataset, which contains tumor microarray data for 128 women with high-grade serous ovarian cancer, where univariate nearest-neighbor G-functions for immune cells have already been extracted. This data has been modified, but the original data and experiment are described in Steinhart et al. (2021). Each subject has one image.
data("ovarian_FDA")
ovarian_FDA
#> mxFDA Object:
#> Subjects: 128
#> Samples: 128
#> Has spatial data
#> Univariate Summaries: Gest
#> Bivariate Summaries: None
#> Multivariate Summaries: None
#> FPCs not yet calculated
#> MFPCs not yet calculated
#> FCMs not yet calculated
#> MFCMs not yet calculated
#> Scalar on Functional Regression not calculatedThese functions are visualized below.
plot(ovarian_FDA, y = "fundiff", what = "uni g", sampleID = "patient_id") +
geom_hline(yintercept = 0, color = "red", linetype = 2) +
theme_minimal() 
Values above the dotted red line are interpreted as having a higher probability of observing a neighboring immune cell at radius \(r\) than would be expected under complete spatial randomness (CSR). Conversely, values below the dotted line exhibit lower probability of observing a neighboring immune cell than expected under CSR.
Run and visualize FPCA
Below we run FPCA on these G-function curves, using the
run_fpca() function, and store these results in the
mxFDAobject called ovarian_FDA. The argument
metric is used to specify which spatial summary function to
estimate FPCA for - in this case "uni g" specifies the
univariate G function, r denotes which variable specifies
the radius, and value denotes which value are the values of
the spatial summary function. The argument pve takes a
number from 0 to 1, and specifies the percent variance explained by the
FPCA decomposition. For example, if pve = 0.98 is selected,
the resulting FPCA will return the number of FPCs that explain 98% of
variance in the data.
ovarian_FDA <- run_fpca(ovarian_FDA,
metric = "uni g",
r = "r",
value = "fundiff",
pve = .95)
#> 128 sample have >= 4 values for FPCA; removing 0 samplesNote that the summary of this object now shows the number of functional principal components (FPCs) that have been calculated:
summary(ovarian_FDA)
#> mxFDA Object:
#> Subjects: 128
#> Samples: 128
#> Has spatial data
#> Univariate Summaries: Gest
#> Bivariate Summaries: None
#> Multivariate Summaries: None
#> FPCs Calculated:
#> Gest: 4 FPCs describe 96.9% variance
#> MFPCs not yet calculated
#> FCMs not yet calculated
#> MFCMs not yet calculated
#> Scalar on Functional Regression not calculatedThe FPCA results are stored in the functional_pca slot
of the ovarian_FDA object. To access this slot directly and
view the extracted summary functions, type:
We can visualize the results of FPCA using the plot()
function (see ?plot.mxFDA for details). The output of
plot() is a {ggplot2} object which can then be
easily added to/manipulated as any ggplot plot would. Below we plot the
first two functional principal components from our FPCA analysis of the
G summary functions.
p1 = plot(ovarian_FDA, what = 'uni g fpca', pc_choice = 1)
p2 = plot(ovarian_FDA, what = 'uni g fpca', pc_choice = 2)
ggarrange(p1, p2, nrow = 1, ncol = 2)
The plots above show \(\hat{\mu}(r)\pm \sqrt{\hat{\lambda}_1\hat{\psi_k}(t)}\) for the first and second FPC. The black line in each plot is the mean of the data- this shows that the peak in probability of observing a neighboring immune cell occurs, on average, at about radius \(r = 12\). In these plots, the blue and red lines represent one standard deviation of the FPC added (blue) or subtracted (red) from the population mean. For example in the plot for FPC 1, we can interpret the pattern of this FPC as a shift up or down from the population mean. Subjects who have a high score, \(c_{i1}\) for FPC 1 have, across radii, a larger probability of observing a neighboring immune cell than average in the dataset. This pattern explains 76.2% of variance in the data.
We can further explore the results of FPCA interactively using the
{refund.shiny} package from Wrobel
et al. (2016). Execute the code below to launch a Shiny app with
interactive plots.
