Statistical analysis

All analyses were performed using open-source R software (version 3.1.3).⁴⁶ Cross-sectional demographic differences between the patients and controls were tested using χ² tests for categorical variables (gender, handedness and ethnicity) and F tests or independent t-tests for continuous variables (age, intracranial volume (ICV), CPZ and time between scans). Longitudinal change in clinical variables was assessed using paired t-tests (duration of illness, PANSS scores and CPZ). Differences in age, ICV, CPZ, duration of illnesses between the data sets were assessed using the Welch two-sample t-test of unequal variances.

We first determined whether there were any group-based differences in the overall mean hippocampus volume of Data set 1 and Data set 2, as well as the early course patients of Data set 1 and matched controls. A multiple linear regression, with volume as the dependent variable, diagnosis as the main predictor, and ICV, age and gender as covariates, was conducted. We then investigated whether schizophrenia differentially affects the volume measures of the inter-related hippocampal subfields. The Shapiro–Wilk test and Bartlett's test of homoscedasticity was first performed to ensure multivariate normality of the subfield volumes and equal variances of in the healthy control and patient groups of each data set, respectively. The subfield values were then log transformed. A multivariate analysis of covariance, with the 14 subfields as dependent (and correlated) variables, diagnosis as the main predictor, and ICV, age, gender and duration of illness as covariates, was conducted. The alpha was set at P<0.05. This was followed by a post-hoc univariate analysis of covariance to determine which subfield (dependent variable) contributed to the significant overall effect of illness.

We sought to determine whether there was an interaction effect between diagnosis and time (between baseline and follow-up scans) for each hippocampal subfield. A separate multi-level model was constructed here, which accounted for the unevenly spaced time-points among subjects and the intra-individual variability in initial subfield volumes and their trajectories (Supplementary Methods 3). After model fitting, fixed effects included diagnosis, time, interaction between diagnosis and time, ICV, CPZ, age, and gender. Random effects included individual intercept and slope of time. The change in volume was modeled linearly, as volume trajectories in studies of gray matter in schizophrenia⁴⁷ and the hippocampus in childhood-onset schizophrenia (which followed subjects until their late twenties)⁴⁸ have been shown to be linear. In addition, the annualized rate of change in subfield volume measures was calculated: (Volume_follow-up – Volume_baseline)/(Volume_baseline × time).

As some studies have reported effects of treatment with antipsychotics,¹¹ antidepressants⁴⁹ and mood stabilizers⁵⁰ on hippocampal structure, we repeated our primary analyses with the dosage or use of medication classes of antipsychotics, antidepressants and mood stabilizers as covariates. To minimize the potential confounding effect of ethnicity, we also repeated our analyses with ethnicity included as a covariate.

Cross-sectional analyses: We examined the relationship between clinical measures (PANSS subscales) and absolute volume measures of the subfields that were significantly different in the patients. A linear regression model was used, with the clinical measures as primary variables of interests, and age and sex among the covariates. We also tested the hypothesis that duration of illness is correlated with subfield volumes. In addition, we tested for correlations between antipsychotic dosages and subfield volumes.

Longitudinal analysis: We examined whether there was an intra-individual relationship between the rate of change in clinical symptoms—calculated by (Scores_follow-up–Scores_baseline)/(Scores_baseline × time)—and the annualized rate of change in the subfield volume measures using similar regression modeling, controlling for age, sex, CPZ and baseline duration of illness.

To address the multiple testing for the various hypotheses, the Holm–Bonferroni method controlling for family-wise errors at alpha level (0.05) was applied.⁵¹