Cord Serum Lipidome in Prediction of Islet Autoimmunity and Type 1 Diabetes

The subjects in this study were chosen from the ongoing prospective DIPP study (14) in which infants in three university hospitals in Finland (Turku, Tampere, and Oulu) are screened for T1D-associated HLA genetic risk alleles (15). Families of children recognized to have increased HLA-conferred risk for T1D are invited to join the study (14). According to the DIPP study protocol, the children are prospectively observed at 3- to 12-month intervals until age 15 years or until the development of clinical T1D. Levels of T1D-associated autoantibodies (islet cell antibodies [ICAs], insulin autoantibodies, IA-2 autoantibodies, and GAD autoantibodies) are determined from the serum samples taken at each follow-up visit (16,17).

This study included DIPP children born between 1994 and 2006 in Turku. According to the DIPP data collected by 31 March 2009, the following groups of children with available cord blood samples were included in the analyses: infants who later developed T1D-associated autoantibodies and progressed to T1D during the follow-up (progressors; N = 33), infants who later developed three or four autoantibodies but have remained clinically unaffected (N = 31), infants who later developed two autoantibodies but have remained clinically unaffected (N = 31), infants who later developed one autoantibody during the follow-up but have remained clinically healthy (N = 48), and clinically unaffected healthy control children (N = 143) without any autoantibodies matched for each child in the other study groups by sex, HLA-conferred T1D risk genotype, and date of birth. All control children were persistently autoantibody-negative from birth and for at least 12 months after the age when the matched case child developed clinical T1D or seroconverted to positivity for the last detected autoantibody during the follow-up. In total, cord serum samples from 286 infants were analyzed in this study (Table 1). Serum of all cord blood samples was isolated within 3 h after collection and stored at −75°C until analyzed.

The Ethics Committee of the Hospital District of Southwest Finland approved the study, and written informed consent was obtained from all participating families.

The analysis of the cord blood HLA-DQB1 genotype was performed with time-resolved fluorometry-based assay using four lanthanide chelate–labeled sequence-specific oligonucleotide probes detecting DQB1*02, DQB1*03:01, DQB1*03:02, and DQB1*06:02/3 alleles (18). Children with T1D risk–associated DQB1*02/DQB1*03:02 or DQB1*03:02/x genotypes (with x representing lack of DQB1*02 and protection-associated DQB1*03:01, DQB1*06:02, and DQB1*06:03 alleles) were eligible for the follow-up program.

The ICAs were detected with the use of indirect immunofluorescence, whereas the three biochemical autoantibodies were quantified with the use of specific radiobinding assays (3). We used cutoff limits for positivity of 2.5 Juvenile Diabetes Foundation units for ICA, 3.48 relative units for insulin autoantibodies, 5.36 relative units for GAD autoantibodies, and 0.43 relative units for IA-2 autoantibodies. The disease sensitivity and specificity of the ICA assay were 100 and 98%, respectively, in the fourth round of the international workshops on standardization of the ICA assay. The disease sensitivity and specificity of the insulin autoantibodies assay were 58 and 100%, respectively, in the 2005 Diabetes Antibody Standardization Program (DASP) workshop. The corresponding characteristics of the GAD autoantibodies assay were 82 and 96%, and those of the IA-2 autoantibodies assay were 72 and 100%. Maternal antibodies that were transferred through the placenta were disregarded in the data analysis.

The global lipidomics approach based on ultra-performance liquid chromatography coupled to mass spectrometry was applied (19). A first standard mixture (20 µL) containing lysophosphatidylcholine (LPC)(17:0), phosphatidylcholine (PC)(17:0/17:0), phosphatidylethanolamine (PE)(17:0/17:0), and ceramide (CER)(d18:1/17:0) (Avanti Polar Lipids, Alabaster, AL) and triglyceride (TG)(17:0/17:0/17:0) (Larodan Fine Chemicals AB, Malmö, Sweden) was added to 10 µL of serum samples, and the samples were extracted with chloroform and methanol (2:1; 100 µL). The lower phase (60 µL) was collected after centrifugation and 20 µL of the second internal standard mixture [LPC(16:1D₃), PC(16:1/16:1-D₆), and TG(16:0/16:0/16:0-¹³C3)] was added. The extracts (2.0 μL) were analyzed on a Waters Q-Tof Premier mass spectrometer (Waters, Milford, MA) combined with an Acquity UPLC (Waters) using an Acquity UPLC BEH C18 2.1- × 100-mm column with 1.7-μm particles. The solvent system included 1% 1 mol/L NH₄Ac, 0.1% HCOOH and acetonitrile/isopropanol (1:1, 1% 1 mol/L NH₄Ac, 0.1% HCOOH) in gradient elution mode with a flow rate of 0.4 mL/min. The lipid profiling was performed using electrospray ionization in positive mode at a mass range of charge/mass ratio 300–1,200 with scan duration of 0.2 s. The data processing using MZmine 2 (20) included detection and alignment of peaks, peak integration, normalization, and peak identification. Only the lipids detected in ≥66% of the samples after the MZmine 2 peak detection were retained for the data analysis. Lipids were identified using an internal spectral library or with tandem mass spectrometry (19). The data were normalized using one or more internal standards representative of each class of lipid present in the samples as previously described (19).

The established lipidomics platform applied in this study covers the major serum lipid classes such as PCs, sphingomyelins (SMs), CERs, PEs, and TGs. However, because the mass spectrometry analysis is performed with electrospray ionization in positive ion mode, the chosen approach is not optimal for the detection of less abundant negatively charged serum lipids such as phosphatidylserines, phosphatidylglycerols, and phosphatidylinositols.

Data were analyzed using R software (21). Lipid concentration values were first log₂-transformed. Clustering analyses were performed using the MCLUST (22) method. In univariate analyses, the log fold changes for every pair of groups were computed as the differences of the mean intensity values between each case and the pooled control group. The 95% CI was computed based on the average SEs of intensity values in each group. Fold changes were obtained by transforming back the log fold changes. Two-sample, two-sided t tests and P value calculations were performed using functions t.test and p.adjust, respectively. The odds ratio (OR) for each variable was computed using the logistic regression. Pearson correlation coefficients and P values were computed for every pair among selected clinical variables and lipid clusters (LCs) using the cor.test function. The results were visualized using the package corrplot.

From each LC, the top 10% of the lipid variables were selected based on the lowest P value. A logistic regression model was built on selected variables using function glm. To find the best submodel with fewer variables, a stepwise backward–forward procedure with Akaike information criterion was applied using function stepAIC from the package MASS (23). To test whether the model was well-calibrated, the Hosmer-Lemeshow goodness-of-fit test was applied (24). To estimate the generalization performance of the model obtained, a cross-validation procedure with 1,000 iterations was applied as follows. In each iteration, a subset containing two-thirds of the samples was randomly drawn for use as training data in such a way that a case and its matching control were always drawn together. The remaining one-third of samples served as the test data. A logistic regression model with the selected variables was trained using the training data and applied to predict the test data. Using the predicted value for each subject on the test data, a receiver operating characteristic curve was constructed using the packages ROCR (25) and pROC (26) to describe classification performance of the logistic regression model. Then, the area under the receiver operating characteristic curve (AUC) was computed. The AUC can be interpreted as rescaled statistic of the Wilcoxon-Mann-Whitney test (27). The optimal cutoff point on the receiver operating characteristic curve was determined to maximize the harmonic mean of sensitivity and specificity pairs. After 1,000 iterations, the mean values of AUC and OR were computed along with their CIs.

The umbilical cord serum lipidome was analyzed from infants who later developed T1D (N = 33), infants who developed three or four islet autoantibodies (N = 31), two islet autoantibodies (N = 31), or one islet autoantibody (N = 48) during the follow-up and controls (N = 143) matched for sex, HLA-DQB1 genotype, city of birth, and period of birth (Table 1). A total of 159 molecular lipids were measured across all 286 serum samples using the established ultra-performance liquid chromatography coupled to mass spectrometry platform (19), of which 137 lipids were identified. We first applied clustering to decompose the lipidomics data into 10 LCs (Table 2). As expected, the division of clusters to a large degree follows different lipid functional or structural groups. For example, the most abundant serum phospholipids including PCs and SMs were included in LC2 and LPCs were included in LC3, whereas the TGs were in LCs 8–10, with the cluster assignments depending on the fatty acid composition of TG species (Table 2).

The lipidomics data were first analyzed at the cluster level (Fig. 1 and Table 2). Pairwise correlation analysis showed that among the 10 LCs, only LC2 was significantly inversely correlated with the study groups (variable study group in Fig. 1A; i.e., low lipid concentrations associate with increased risk of T1D). This association is specific for T1D progressors, whereas no significant changes were observed in other study groups (Fig. 1B). None of the LCs were correlated with HLA-associated T1D risk (variable HLA risk in Fig. 1A). A pattern similar to that of the cluster level was observed at the molecular lipid level. The top 12 ranked lipids based on odds of progression to T1D were from LC2 (Table 3), including the most abundant SMs and PCs, such as SM(d18:1/24:1), SM(d18:1/16:0), PC(16:0/18:1), and PC(18:0/18:1).

Because it has been previously documented that age at seroconversion for islet autoantibody positivity is associated with different metabolic profiles before and after seroconversion (12), we also investigated whether age at seroconversion is associated with the cord serum lipidomic profiles. When age at seroconversion was treated as a continuous variable, no significant associations were found for any of the LCs (Fig. 1A). When dividing the T1D progressors into two groups by age at seroconversion (early autoantibody-positivity indicated by seroconversion at age younger than 4 years and late autoantibody positivity indicated by seroconversion at age older than 4 years), LC2 lipids were significantly decreased only in the early autoantibody-positive group (Fig. 1C). However, one cannot rule out that this is attributable to the lack of statistical power because the late autoantibody-positive group included only 10 T1D progressors and the average degree of decreased LC2 was similar to that in the early autoantibody-positive group (Fig. 1C).

Among other clinical variables, gestational age was positively associated with several LCs including all three TG clusters (Fig. 1A). At the individual lipid level, the strongest associations were observed for TGs from LC9 and LC10, i.e., TGs that contain essential fatty acids. Out of 30 top-ranked lipids according to the P value (null hypothesis [H₀: R = 1]), 19 were from LC9 and 6 were from LC10 (P < 10⁻⁸ for all of them). None of these top-ranked lipids were significantly associated with age at diagnosis or age at initial seroconversion. These findings are in agreement with previous studies in which the fatty acid composition of TG fractions was assessed in pregnant women and infants (28,29).

To assess the feasibility of predicting T1D, we performed a model selection in multiple cross-validation runs as described. The diagnostic model was developed based on the data from T1D progressors and their matched controls. The best model derived from logistic regression analysis was obtained by combining seven lipid metabolites (Fig. 2A). Tests for goodness-of-fit revealed no evidence for the lack of fit of the model (P = 0.43), indicating that it was well-calibrated. The model had an AUC of 0.71 (95% CI, 0.53–0.85) and an OR of 5.94 (95% CI, 1.07–17.50). When using a different set of controls to test the model (those matched to the three or four autoantibodies group, two autoantibodies group, and one autoantibody group), the OR was 2.71 (95% CI, 1.08–6.79). The prediction of T1D was poor in the three or four autoantibodies group, two autoantibodies group, and one autoantibody group, with the ratio of percent predicted T1D cases in each group compared with percent predicted cases in pooled controls from the three groups being 1.00, 0.84, and 0.96, respectively.

Next, we sought to predict development of β-cell autoimmunity. Among the different case groups, seroconverters for one autoantibody were least likely to develop T1D (16,30) and thus independent from T1D progression. Therefore, we developed a model to predict exclusively β-cell autoimmunity using the data from infants in the one autoantibody group and their matched controls. The best model was obtained by combining seven lipid metabolites (all of the lipids were different from the ones included in the lipid signature to predict T1D; Fig. 2B). Tests for goodness-of-fit revealed no evidence for the lack of fit of the model (P = 0.63), indicating that it was well-calibrated. The model had an AUC of 0.70 (95% CI 0.57–0.83) and an OR of 4.54 (95% CI 1.07–12.00). When using a different set of controls to test the model (those matched to the T1D, three or four autoantibody, and two autoantibody groups), the OR was 5.11 (95% CI 2.40–10.87). The lipid signature associated with progression to one autoantibody was identified in 25% of all controls, in 68% of infants in the one autoantibody group, in 55% in the two autoantibody group, in 48% in the three or four autoantibody group, and in 39% in the T1D group.

For both predictive models, inclusion of other clinical variables such as birth weight and gestational age did not improve the models: ORs were 4.86 (95% CI 0.8–12.8) for the T1D model and 4.07 (95% CI 1.07–10) for the model to predict β-cell autoimmunity.

Our study confirms previous findings (11) that progression to T1D is associated with decreased concentrations of major choline-containing phospholipids (SMs and PCs) in cord blood. Here, we also show that the reduction in phospholipids is specifically associated with progression to T1D but not with β-cell autoimmunity in general and is predominantly seen in T1D children who seroconverted before age 4 years.

Low PC and SM levels in T1D progressors cannot be explained by genetics alone because the subjects in our study were matched with the controls for the single most important gene locus, HLA-DQB1. Choline is an essential nutrient that is mainly found in plasma and not in free-form as part of the head groups of SMs and PCs. Choline is in particularly high demand during pregnancy as a substrate for building cellular membranes because of rapid fetal tissue expansion and for increased production of lipoproteins (31). Furthermore, choline is a major donor of methyl groups needed for DNA methylation and therefore is essential for developmental processes, including genomic imprinting and the maintenance of genome stability (31,32). For example, maternal choline intake during pregnancy alters the epigenetic state of fetal cortisol-regulating genes (33). In a recent study of pregnant women, it has been shown that increased dietary choline intake during pregnancy increases its utilization as a methyl donor in both maternal and fetal compartments (34). Although dietary supplementation mainly has been considered as a way to increase the levels of choline, other factors affecting the choline metabolism such as gut microbiota (35) also may need to be considered to optimize the utilization of dietary choline. Together, whereas the exact causes for the observed metabolic profile in T1D progressors at birth remain unknown, our findings do suggest an intriguing possibility that environmental factors during pregnancy affect the choline metabolism and, hence, genomic imprinting of the fetus, leading to increased T1D risk.

Our findings may have wider clinical significance beyond T1D because they suggest that lipidomics, or more broadly metabolomics, may be useful in newborn screening not only for the risk of in-born errors of metabolism (36) but also for the risk of multifactorial chronic diseases. Lipidomics technology has been maturing rapidly (13,37), and translation of findings such as ours in the current study into an inexpensive assay applicable in the health care setting is already feasible.

In conclusion, our study suggests that decreased choline-containing phospholipids in umbilical cord blood are predictors of T1D but not of β-cell autoimmunity per se. The lipidomic profiles associated with progression to T1D may have diagnostic value and may help identify infants at risk for development of T1D. Our findings primarily implicate altered choline metabolism in pregnancy as a contributing factor behind the observed lipidomic profiles in cord sera of T1D progressors. Future studies focusing on metabolism and genomic imprinting during pregnancy and neonatal life will be important to disclose the character of gestational events contributing to reduced phospholipid levels in relation to the development of T1D in the offspring.

The study was supported by Academy of Finland (Centre of Excellence in Molecular Systems Immunology and Physiology Research 2012-2017, Decision No. 250114, to M.O., M.K., and O.S.), the DIAPREPP project of the Seventh Framework Program of the European Community (HEALTH-F2-2008-202013 to M.O. and O.S.), and Juvenile Diabetes Research Foundation (36-2008-919 to M.O. and O.S.).

No potential conflicts of interest relevant to this article were reported.

M.O. was responsible for the study concept and design, analyzed and interpreted data, drafted the manuscript, critically revised the manuscript for important intellectual content, and supervised the study. P.G. and J.M. analyzed and interpreted data. N.L., M.M., H.N., S.S., and V.S. were responsible for acquisition of data. H.H. and R.V. were responsible for critical revision of the manuscript for important intellectual content. J.I. and M.K. were responsible for acquisition of data and critical revision of the manuscript for important intellectual content. M.S.-A. analyzed and interpreted data, was responsible for acquisition of data and critical revision of the manuscript for important intellectual content, and supervised the study. T.H. was responsible for acquisition of data and critical revision of the manuscript for important intellectual content and supervised the study. O.S. was responsible for the study concept and design and critical revision of the manuscript for important intellectual content and supervised the study. M.O. is the guarantor of this work and, as such, had full access to all of the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.

The authors thank the dedicated personnel of the DIPP study in Turku, Oulu, and Tampere and the DIPP families for their support. The authors acknowledge the technical support in lipidomics analysis from Ulla Lahtinen, Anna-Liisa Ruskeepää, Sandra Castillo, Erno Lindfors, and Tuulikki Seppänen-Laakso (VTT Technical Research Centre of Finland, Espoo, Finland).