Longitudinal functional principal component analysis

Abstract

We introduce models for the analysis of functional data observed at multiple time points. The dynamic behavior of functional data is decomposed into a time-dependent population average, baseline (or static) subject-specific variability, longitudinal (or dynamic) subject-specific variability, subject-visit-specific variability and measurement error. The model can be viewed as the functional analog of the classical longitudinal mixed effects model where random effects are replaced by random processes. Methods have wide applicability and are computationally feasible for moderate and large data sets. Computational feasibility is assured by using principal component bases for the functional processes. The methodology is motivated by and applied to a diffusion tensor imaging (DTI) study designed to analyze differences and changes in brain connectivity in healthy volunteers and multiple sclerosis (MS) patients. An R implementation is provided.87.

Fig 1

Top: Sagittal image of the corpus callosum in one of the study subjects, a healthy 33-year-old man, showing the segmentation used [following 46] for construction of the tract profile. Values denote the bin number at the boundary point from the splenium (back of the head) to the genu/rostrum (closer to the eyes). Bottom: Two example subjects (both MS patients) from the tractography data with 5 and 6 complete visits, respectively. Shown are the fractional anisotropy along the corpus callosum, measured at the 120 sample points. Different visits for the same subject are indicated by color and overlaid.

Fig 2

True and estimated eigenfunctions ${\phi }_k^X=({\phi }_k^0,{\phi }_k^1)$ and ${\phi }_k^U$, k = 1, …, 4. The left column gives results for the part ${\phi }_k^0$ corresponding to the random functional intercept X_i,0, the middle column for the part ${\phi }_k^1$ corresponding to the random functional slope X_i,1, and the right column for the component ${\phi }_k^U$ corresponding to the visit-specific functional deviation U_ij. Shown are the true function (thick black line), the mean of the estimated functions over 1000 simulations (dashed red line), the pointwise 5th and 95th percentiles of the estimated functions (blue), and the estimated functions from the first 50 simulations (grey). Simulations were based on model (2.3) with N_X = N_U = 4, a balanced design with I = 100 and J_i = 4 for all i, non-orthogonal ${\phi }_k^0,{\phi }_k^1$ and ${\phi }_k^U$ with unequal weight on ${\phi }_k^1$ and ${\phi }_k^0$, a mixture distribution for the scores ξ_ij and ζ_ijk, and no bivariate smoothing of the covariance functions.

Fig 3

Computation time for LFPCA for a simulated data set with the given number of subjects I and number of observations per subject J, and with D sample points per curve. Specifics of how computation time was measured are given in Section 4.3.

Fig 4

The first three estimated principal components for the random intercept and slope process X. The left column gives estimates for the ${\phi }_k^0$, corresponding to the random functional intercept X_i,0. Depicted are estimates for the overall mean η(d) (solid line), and for $\eta (d)\pm 2\sqrt{{\lambda }_k}{\phi }_k^0$, k = 1, 2, 3 (+ and −, respectively). The middle column gives the corresponding results for the random functional slope X_i,1. The right column shows boxplots for the estimates of the scores ξ_ik corresponding to ( ${\phi }_k^0,{\phi }_k^1$), k = 1, 2, 3, by case/control group. Estimated scores for the two example patients with tract profiles shown in Figure 1 are indicated by A and B, respectively.

Fig 5

The first three estimated principal components for the visit-specific deviation process U. The left column gives results for the principal components ${\phi }_k^U$, depicting estimates for the overall mean η(d) (solid line), and for $\eta (d)\pm 2\sqrt{{\nu }_k}{\phi }_k^U$, k = 1, 2, 3 (+ and −, respectively). The right column shows boxplots for the estimates of the scores ζ_ik corresponding to ${\phi }_k^U$, k = 1, 2, 3, by case/control group. Estimated scores for example visits of the two patients with tract profiles shown in Figure 1 are indicated by A (visit 8) and B (visit 2), respectively.