Samples

Structural brain images from 881 subjects were collected with T1-weighted magnetic resonance imaging (MRI) from three independent cohorts: a sample collected at San Antonio (the SA sample), a sample collected at Nathan Kline Institute - Rockland (the NKI sample)66 and a sample collected by the Centers of Biomedical Research Excellence (the COBRE sample). The latter two samples were publicly available. All patients completed the Structured Clinical Interview for DSM-IV Axis I Disorders (SCID) and met the corresponding DSM-IV criteria for MDD, bipolar 1 disorder (BD-I) or SCZ. Healthy controls (HC) had no current condition or history of psychiatric dysfunctions. All subjects had no history of neurological diseases. Signed informed consent have been obtained from all the subjects. No identifiable image from a specific subject was used for publication. The study protocols were approved by the Institutional Review Boards at the University of Texas at San Antonio in accordance with their guidelines and regulations.

Image Acquisition and Processing

The SA sample was all acquired at one Philips 1.5 Tesla MRI scanner (Philips Medical System, Andover, MA, USA) with a three-dimensional axial fast field echo sequence with the following parameters: repetition time (TR) = 24 ms, echo time (TE) = 5 ms, flip angle = 40°, field of view (FOV) = 256 mm, slice thickness = 1 mm, matrix size = 256 × 256 and 150 slices. All scans were visually inspected to rule out gross artifacts. The NKI sample was acquired at one 3 Tesla Siemens Magnetom TrioTim syngo scanner with a three-dimensional magnetization-prepared radio-frequency pulses and rapid gradient-echo (MPRAGE) sequence with the following parameters: TR = 2500 ms, TE = 3.5 ms, flip angle = 8°, FOV = 256 mm, slice thickness = 1 mm, matrix size = 256 × 256 and 192 slices (see http://fcon_1000.projects.nitrc.org/indi/pro/nki.html). The COBRE sample was acquired at one 3 Tesla Siemens Magnetom TrioTim syngo scanner with a three-dimensional multi-echo (ME-MPRAGE) sequence with the following parameters: TR = 2530 ms, TE = 1.64, 3.5, 5.36, 7.22, 9.08 ms, flip angle = 7°, FOV = 256 mm, slice thickness = 1 mm, matrix size = 256 × 256 and 176 slices (see http://fcon_1000.projects.nitrc.org/indi/retro/cobre.html).

The structural brain images were preprocessed and the cortical surface of each brain was reconstructed with Freesurfer67, 68 (version 5.3, http://surfer.nmr.mgh.harvard.edu/). The reconstructed images were visually inspected by the authors to exclude the apparent reconstruction errors. The local gyrification index (GI) at each vertex of the reconstructed cortical surface mesh was calculated using the toolbox in Freesurfer with default settings69. Briefly, a circular region of interest was delineated on the outer surface, and its corresponding region of interest on the inner cortical surface was identified using a matching algorithm as described elsewhere70 and the ratio between the folded inner cortical surface and its corresponding exposed outer surface was calculated as GI. The resulted GI values of each subject were smoothly sampled (Gaussian kernel, 10 mm) onto an average template provided by Freesurfer (fsaverage; 163842 vertices), so that cross individual comparison could be performed. The averaged GI for each subject was also calculated to represent the gyrification of whole brain.

The effect of brain volume on GI was taken into account for by adjusting the GI values with a linear regression of intracranial volume (ICV) on GI at each vertex and the whole brain GI, because previous studies found that the GI could be correlated with the brain volume31, 71.

Modeling Gyrification Trajectory as a function of Age

We compared the fitting of the whole brain GI with age from ten common mathematical functions and chose the best function based on the averaged mean squared errors (MSE). The MSE of each fitting function was calculated with the non-linear fitting functions of the Statistics and Machine Learning Toolbox in MATLAB (The MathWorks, Inc., Natick, Massachusetts, United States) using default settings. Because the COBRE sample only included adults, it was not involved in the function selection process. For SA and NKI samples, the MSE of each function was calculated based on the residuals derived from two types of validation procedures: (a) a leave-one-out cross-validation within each of the two samples and (b) a cross-sample validation, in which a fitting function was optimized for one sample and tested on the other. Thus, we had four MSEs for each function and the averaged of these MSEs was used as the fitting performance for the corresponding function as shown in Table S1. The best function was the logarithmic function with three parameters.

As a result, we used the following logarithmic function of age to fit the average GI trajectory:

$${\rm{GI}}=a+b\ast ln(age+c)$$ (1)

where a was an indicator of GI levels independent of age, b controlled the decrease rate of GI over age and c was the translational term of age. The initial value of a was set to 4, because the mean of the brain average GI values of all the HC subjects was around 3 and the maximum whole brain GI was near 4. The initial value of b was set to −1, because GI apparently decreased with age. The initial value of c was set to 0. We found that c was sensitive to the age range of the sample and we had varied age range in COBRE HC, MDD, BD and SCZ samples. Furthermore, c did not provide any information of the GI levels and trajectory shapes, and its value was consistent in SA and NKI HC samples. In order to reliably compare the GI trajectories of all the samples without changing the shape and levels of GI trajectory, we fixed the c value to −2.9991 according to the cross-validated fitting for the combined SA and NKI HC sample. Then each of the MDD, BD and SCZ groups was fitted separately.

GI was adjusted for study site by an amount of the estimated GI difference between the HC in each cohort and the combined SA and NKI HC sample at a given age [a HC + b HC \(\ast \) ln(age + c)] − [a X + b X \(\ast \) ln(age + c)], where X represents the sample of HC in SA, NKI and COBRE cohorts and the mean adjustment values were as small as −0.0023, 0.0043 and 0.0523, respectively. The adjustment for the COBRE cohort was larger than the other two cohorts, because the GI in the COBRE HC sample was generally lower than the HC in the other two cohorts, which could also be observed in the mean GI of adult HC in the COBRE cohort (2.748) compared to the GI of adult HC in SA (2.7938) and NKI cohorts (2.7777). To further confirm whether the adjustment introduced significant effect in the SCZ sample in the COBRE cohort, we also compared the GI over the brain of HC and SCZ within the COBRE cohorts.

Although the samples showed different distributions of males and females (see also72), which might have an effect on the GI trajectory of the samples with major psychiatric disorders, we found that effect of gender on GI was negligible after GI was adjusted by the ICV (Supplementary Materials Fig. S1). This might indicate that the gender effect on GI could be mostly explained by the brain volume. Thus, no further adjustment for gender was performed.

Estimating the Parameters with Resampling

In order to quantify the difference between the GI trajectories of HC, MDD, BD-I and SCZ, it was necessary to estimate the distributions of fitting parameters. We utilized the resampling technique. Briefly, for each of the HC, MDD, BD and SCZ samples, the sample was re-sampled 10000 times and the fitting was performed for each of the 10000 samples. In order to stratify the age distribution during the resampling, each sample was divided into four age blocks: 4–9, 9–18, 18–40 and 40–83. During each iteration of the resampling, 50% of the subjects in each age block were selected without duplicate. Thus, we had 10000 sets of parameters for each sample, and it was then possible to estimate the distributions of the fitting parameters in different samples.