Polygenic Multiple Sclerosis Risk and Population-Based Childhood Brain Imaging.

OBJECTIVE: Multiple sclerosis (MS) is a neurological disease with a substantial genetic component and immune-mediated neurodegeneration. Patients with MS show structural brain differences relative to individuals without MS, including smaller regional volumes and alterations in white matter (WM) microstructure. Whether genetic risk for MS is associated with brain structure during early neurodevelopment remains unclear. In this study, we explore the association between MS polygenic risk scores (PRS) and brain imaging outcomes from a large, population-based pediatric sample to gain insight into the underlying neurobiology of MS.
METHODS: We included 8- to 12-year-old genotyped participants from the Generation R Study in whom T1-weighted volumetric (n = 1,136) and/or diffusion tensor imaging (n = 1,088) had been collected. PRS for MS were calculated based on a large genome-wide association study of MS (n = 41,505) and were regressed on regional volumes, global and tract-specific fractional anisotropy (FA), and global mean diffusivity using linear regression.
RESULTS: No associations were observed for the regional volumes. We observed a positive association between the MS PRS and global FA (β = 0.098, standard error [SE] = 0.030, p = 1.08 × 10-3 ). Tract-specific analyses showed higher FA and lower radial diffusivity in several tracts. We replicated our findings in an independent sample of children (n = 186) who were scanned in an earlier phase (global FA; β = 0.189, SE = 0.072, p = 9.40 × 10-3 ).
INTERPRETATION: This is the first study to show that greater genetic predisposition for MS is associated with higher global brain WM FA at an early age in the general population. Our results suggest a preadolescent time window within neurodevelopment in which MS risk variants act upon the brain. ANN NEUROL 2020;87:774-787.

Multiple sclerosis (MS) is a severe neurological disorder caused by demyelination in the central nervous system (CNS).1 In addition to distinctive demyelinating lesions, several extensive brain differences have been reported in adult and pediatric MS patients. For example, brain imaging studies have shown atrophy of gray matter (GM), particularly in the thalamus, hippocampus, and other subcortical structures,2, 3 as well as alterations in the microstructure of normal‐appearing white matter (WM).2, 4 Studies in children with a diagnosis of MS have shown reduced total and regional brain volumes at first presentation of symptoms compared to their healthy peers,5 which raises the question of whether neurodevelopmental differences are present before the clinical presentation of MS.

Genetic factors have been shown to contribute considerably to MS pathophysiology, as is evidenced by an increased relative risk (9.5–116.7 fold) of MS in individuals with an affected first‐degree family member and heritability estimates of approximately 30%.6, 7 In addition, previous work has shown evidence of WM lesions in asymptomatic relatives of adult MS patients.8 A recent genome‐wide association study (GWAS) has identified 233 genome‐wide significant genetic variants (single nucleotide polymorphisms [SNPs]), mainly with immunological functions, associated with the risk of MS.9 This work demonstrated that the risk of MS is highly polygenic, with a large number of common genetic variants each contributing a small effect to the overall MS disease risk.

One technique used to quantify the effect of combined genetic risk of multiple common risk variants across the genome is the use of polygenic risk scores (PRS).10, 11 An increased PRS of MS has shown to be associated with a higher risk of MS and an increased susceptibility to pediatric‐onset MS.12, 13 In addition, PRS can be used to explore whether a high genetic risk leads to a higher degree of disease‐related alterations, such as in brain imaging measures. Large population‐based imaging studies in adults from the general population show that smaller GM volume is related to higher genetic risk for MS.14 Further, higher fractional anisotropy (FA) and lower mean diffusivity (MD) were observed in several WM tracts for a PRS for MS, as well as numerous individual genetic risk variants, although not passing statistical significance testing correction.14, 15

Importantly, most previous studies investigating the PRS for MS focused only on those variants that passed genome‐wide significance (except for Brown et al15), omitting other common risk variants that likely also contribute substantially to the risk of MS.9 Furthermore, associations with brain imaging outcomes have not been studied in younger populations, which may provide valuable insight into the link between genetic risk of MS and neurodevelopment during childhood.

Within this context, it is our goal to study the association between genetic risk for MS using PRS and brain magnetic resonance imaging (MRI) that indexes neurodevelopment in school‐age children, including volumetric and microstructural brain measures. As the window of opportunity for prevention and treatment of MS may exist during childhood, we aimed to assess whether genetic risk may index brain alterations earlier than has been expected.

Subjects and Methods

Study Sample

This study was conducted using data from the Generation R Study, a longitudinal population‐based birth cohort (n = 9,749) designed to examine different facets of pediatric development.16 Between March 2013 and November 2015, 3,992 participants aged 8 to 12 years visited a study‐dedicated research center for an MRI scan.17 The present study included unrelated participants of European ancestry with usable MRI and genotype data (Fig 1). Children with perinatal WM abnormalities, as assessed by a neuroradiologist, were excluded. None of the included participants showed other WM abnormalities such as radiologically isolated syndrome.

Figure 1

Flowchart describing the selection process of the study. DTI = diffusion tensor imaging; QC = quality control; T1‐W = T1‐weighted.

To replicate our findings, we identified a sample of participants included from an earlier phase of the Generation R imaging study who did not participate in the current study (n = 186; for study details, see White et al18 and Muetzel et al19). The Generation R Study protocol was approved by the Medical Ethical Committee of the Erasmus Medical Center; the study was conducted according to the Declaration of Helsinki. Written informed consent was obtained from the children's parents or legal representatives.

Genotype Data

Genotype data were obtained from cord blood at birth or from venipuncture during a visit to the research center using either a 610K or 660K SNP array (Illumina, San Diego, CA). Sample collection and genotype calling procedures have been described elsewhere.20 Information about additional quality control procedures of the genotype data, genotype imputation, and calculation of principal components has been reported in a previous study.21 In short, we included participants of European ancestry based upon the first 4 principal components inside the range of the HapMap Phase II Northwestern European founder population and used the 1,000 Genomes (phase I version 3) imputed genotype data for the calculation of our PRS.22, 23

Polygenic Scores

We used a recent large discovery GWAS of MS in 41,505 individuals (14,802 cases and 26,703 controls)9 to compute a weighted PRS based on imputed genotype data. Results from the MS discovery GWAS were obtained via the International Multiple Sclerosis Genetics Consortium (IMSGC; http://imsgc.net/publications/). PRS was calculated using PRSice 2,24 an R‐script that carries out PRS calculation across several p value thresholds (PTs) using PLINK (v1.9).25, 26

We removed high linkage disequilibrium regions (r 2 < 0.10, kB window = 250kB) by using the Generation R European sample as a reference. After this, we used several PTs for inclusion of SNPs in the computation of the PRS (PT < 0.001, 0.005, 0.01, 0.05, 0.1, 0.5, 1). We also calculated a PRS including only SNPs that reached genome‐wide significance (PT < 5 × 10−8) in the final meta‐analysis from the IMSGC (n = 115,803; 47,429 cases and 68,374 controls).9

Given the strong effect of variants located in the major histocompatibility complex (MHC), we calculated 2 sets of PRS, 1 including all available SNPs and another excluding SNPs from the MHC region. For characterization of HLA‐DRB1*15:01, rs3135388 was used as the tag SNP.27

The number of SNPs included in a given PRS at different PTs is shown in Table 1. The distribution of the calculated PRS across the whole study population can be found in Figure 2. The PRS were comparable across gender and independent of age at scan (see Fig 2).

Table 1

Total Number of SNPs in the Different PRS

PRS	PT < 5 × 10−8	PT < 0.001	PT < 0.005	PT < 0.01	PT < 0.05	PT < 0.1	PT < 0.5	PT < 1
All available variants, n	198	1,854	4,963	7,920	24,593	40,615	118,343	161,270
Excluding MHC variants, n	191	1,704	4,800	7,751	24,399	40,406	118,101	161,012

Amount of single nucleotide polymorphisms (SNPs) included in the polygenic risk score (PRS) at each p value threshold (PT).

MHC = major histocompatibility complex.

Figure 2

Descriptive characteristics of the polygenic risk scores at the different p value thresholds (PTs) across the whole study population (n = 1,259). (A) Density distribution of the polygenic risk scores across the whole study population at the different PTs. (B) Distribution of the polygenic risk scores at the different PTs across gender. (C) Height of the polygenic risk scores across the age span of our study at the different PTs.

Brain

Imaging of the brain was performed on a single study‐dedicated 3T MR750w Discovery scanner (General Electric, Milwaukee, WI). Both high‐resolution structural MRI and diffusion‐weighted images were collected. A detailed description of the scan protocol, imaging procedures, and subsequent processing of the imaging data can be found in earlier work from our study group.17, 19, 28

For analyses involving volumetric brain measures, we used T1‐weighted images. Regions of interest (ROI) were defined (including thalamus, hippocampus, total subcortical volume, total GM, and total WM) using the FreeSurfer software package (v6.0)29 and the Desikan–Killiany atlas, which contains 34 cortical parcels per hemisphere.30 For diffusion tensor imaging (DTI) analyses, we used the functional MRI of the Brain's Software Library (FMRIB FSL) and the Camino Diffusion MRI Toolkit to preprocess the data.31, 32 The AutoPtx probabilistic tractography plugin was then used to identify 12 major WM tracts.33 Within these tracts, diffusion characteristics were used to quantify average FA and MD. To model a single latent factor of global FA and MD, representing global indicators of WM microstructure in the brain, these 12 tracts were submitted to confirmatory factor analysis (CFA) using the lavaan R package.19, 28 The WM tracts included in the CFA model, the standardized loading factors, and the fit indices of the global factors for FA and MD resulting from the model can be found in Supplementary Table 1.

The participants in the replication sample were scanned on a 3T MR750 Discovery MRI scanner between September 2009 and February 2012 and were between 6 to 10 years old. During the brain imaging session, both structural and DTI images were acquired. The structural brain data images were segmented using Freesurfer (v5.3),34 and DTI data were processed using the FMRIB FSL and the Camino Diffusion MRI Toolkit.19

Nonverbal

Nonverbal IQ data were studied to assess the association between the MS PRS and cognitive functioning. Nonverbal IQ was assessed at the age of 6 years in around 6,000 participants using 2 subsets of a nonverbal IQ test: Snijders‐Oomen Nonverbal Intelligence Test–Revised (SON‐R 2.5–7; http://www.testresearch.nl). The subsets for Mosaics (spatial visualization abilities) and Categories (abstract reasoning abilities) were used to compute age‐ and sex‐normalized nonverbal IQ scores.

To investigate a possible relationship between the MS PRS, WM microstructure, and cognitive functioning, we performed an additional supplementary analysis between the PRS and nonverbal IQ in participants who had DTI data available.

Statistical Analysis

All analyses were performed using the R statistical software package (v3.5.1).25 Before running each model, all variables of interest were checked for normality. We applied t tests to assess the difference in PRS for MS (PT < 1) and descriptive characteristics (gender and level of maternal education) with the participants that did not participate in the volumetric and/or DTI scanning.

For the volumetric brain analyses, the PRS of MS across different thresholds were regressed on total brain volume (TBV). For DTI, the PRS of MS with different thresholds were regressed on global FA and MD. In secondary analyses, the PRS based on the PT that showed the strongest association with the global outcome in the primary analyses was regressed on the individual regional metrics (ie, ROI volumes for the volumetric analyses and the different WM tracts for DTI analyses). In WM tracts that showed evidence of association with the PRS, we ran post hoc analyses testing the associations between the PRS of MS and radial and axial diffusivity (RD and AD, respectively) to gain further insight into the diffusion profile that contributes to differences in FA or MD.

All analyses were adjusted for age at scan, sex, and the first 10 genetic principal components, the latter to adjust for residual effects due to population substructure. Volumetric analyses for specific brain regions were adjusted for total intracranial volume to mitigate the possible confounding effects of head size. In the DTI analyses, we performed sensitivity analyses in the WM tracts correcting for the global diffusion factor to explore whether the PRS associations are primarily global or tract specific. To correct for multiple testing, we used false discovery rate (FDR)35 on the total number of statistical tests for each risk score, PT, and global, tract‐specific, and regional measures. FDR‐corrected significance threshold was applied separately for DTI and morphology. This translated to an uncorrected PT of 0.013 for the DTI analyses and 0.004 for the volumetric analyses.

Lastly, we performed linear regression analyses using inverse probability‐of‐censoring weights (IPCWs) to account for potential differences between the study sample and the Generation R base sample.36 In short, descriptive characteristics that differed between participants included in the MRI study and those in the Generation R base sample (excluding those who participated in the present study) were fitted as covariates into a logistic regression model in the Generation R base sample to estimate a pseudopopulation with a conditional probability and weight for each subject. These weights were used in the study sample regression analyses to account for possible induced selection bias.37

Results

Sample Selection and Information

T1‐weighted structural MRI scans were completed in 3,850 participants, and DTI scans were available in 3,786 participants. A total of 3,186 participants had eligible T1‐weighted scans after excluding images that were collected using different acquisition parameters (n = 22), those with poor image quality (n = 620), and images including major incidental findings (n = 22). Scanner software was upgraded during the study period, which did not affect T1‐weighted volumetric data. However, DTI data were influenced and participants were excluded based on scanner software version. A total of 3,050 participants were found to have eligible DTI scans after excluding participants with poor image quality (n = 447), different version of scanner software (n = 276), and incidental findings (n = 13). Genotype data were available in 2,040 of the participants with T1‐weighted scans and in 1,903 of the participants with eligible DTI data. Subsequent selection on European ancestry, relatedness, and genotype quality resulted in 1,136 participants eligible for volumetric analyses and 1,088 participants eligible for DTI analyses (see Fig 1).

Participants included in the volumetric analyses had a median age of 9.95 years (interquartile range [IQR] = 9.76–10.31), and sex was evenly distributed (49.4% female; Table 2). The distribution of the different T1‐weighted volumes between gender and across the age span is shown in Figure 3. No difference was found in mean PRS of MS (p = 0.54) or in gender distribution (p = 0.57) compared to participants who did not participate in the volumetric MRI scanning (n = 1,694). However, we observed higher levels of maternal education in the participants eligible for the T1‐weighted MRI analyses (p < 0.01). The median age in the participants included for the DTI analyses was 9.96 years (IQR = 9.78–10.36), with an even distribution of sex (49.2% female; see Table 2). Figure 4 shows the distribution of global FA and MD between the genders and across the different ages. There were no differences in mean PRS for MS (p = 0.96) and gender (p = 0.66) compared to genotyped participants of European ancestry who did not participate in the DTI scanning (n = 1,742). In the sample eligible for the DTI analyses, we again had a relative overrepresentation of high maternal education (p < 0.01).

Table 2

Demographic and Educational Information about the Study Samples

	T1‐W Volumetric MRI, n = 1,136	DTI, n = 1,088	p
Female, n (%)	561 (49.4)	536 (49.2)	0.99
Age at MRI scan, median (IQR), yr	9.95 (9.76–10.31)	9.96 (9.78–10.36)	0.49
Level of maternal education, n (%)
Low‐middle	348 (30.6)	328 (30.1)	0.84
High	767 (67.5)	739 (67.9)	0.87
Unknown	21 (1.8)	21 (1.9)	1.00
Presence of maternal MS, n (%)	6/964 (0.6)	4/920 (0.4)	0.81
Presence of paternal MS, n (%)	2/779 (0.3)	4/758 (0.5)	0.90

DTI = diffusion tensor imaging; IQR = interquartile range; MRI = magnetic resonance imaging; MS = multiple sclerosis; T1‐W = T1‐weighted.

Figure 3

Distribution of the T1‐weighted volumetric outcomes in our study population (n = 1,136). (A) Distribution of the volumetric outcomes across gender (in mm3). (B) Size of our volumetric outcomes across the age span of our study population (in mm3).

Figure 4

Distribution of global fractional anisotropy (FA) and mean diffusivity (MD) in our study population (n = 1,088). (A) Distribution of our diffusion outcomes across gender. (B) Development of our diffusion outcomes across the age span.

Figure 5 shows the PRS–TBV associations across different thresholds and the explained variance in TBV by the PRS. Several different thresholds showed positive associations with TBV (β = 0.055, standard error [SE] = 0.025, ΔR 2= 0.0030, p = 0.03). However, this was not significant after multiple testing correction. Supplementary Table 2 shows the full regression results of the PRS on TBV.

Figure 5

Explained variance in total brain volume (TBV) by the different polygenic risk scores (PRS). (A) The y‐axis represents the increase in explained variance (R 2) by including the PRS of multiple sclerosis in the model. The shades of the bars represent the different p value thresholds (PTs) for inclusion of single nucleotide polymorphisms (SNPs) in the PRS. (B) The y‐axis represents the standardized regression coefficients of the different PRS on TBV, corrected for age, sex, and 10 genetic principal components. [Color figure can be viewed at http://www.annalsofneurology.org]

The PRS PT with the strongest association with TBV was regressed on different subcortical volumes to test for specific regional associations (PT < 0.01), but none were found to be significant. Supplementary Table 3 shows the regression results of the PRS for MS on the different volumes.

The PRS for MS showed significant positive associations with global FA across different PT (Fig 6). The PRS with a PT < 0.01 had the strongest association with global FA (β = 0.098, SE = 0.030, ΔR 2= 0.0095, p = 1.08 × 10−3; Supplementary Table 4). Visual representations of the significant associations with global FA are shown in Figure 6.

Figure 6

Explained variance in global fractional anisotropy (FA) and global mean diffusivity (MD) by the different polygenic risk scores (PRS). (A) The y‐axis represents the increase in explained variance (R 2) by including the PRS of MS in the model. The shades of the bars represent the different p value thresholds (PTs) for inclusion of SNPs in the PRS. (B) The y‐axis represents the standardized regression coefficients of the different PRS on global FA and global MD, corrected for age, sex, and 10 genetic principal components. (C) Scatter plots of the significant associations between our MS PRS and global FA after correcting for multiple testing. [Color figure can be viewed at http://www.annalsofneurology.org]

The PRS showed mainly negative associations with MD (β = −0.050, SE = 0.029, ΔR 2= 0.0025, p = 0.09), but these were not significant (see Fig 6, Supplementary Table 5). The analyses using PRS of MS excluding the HLA‐DRB1*15:01 tag variant rs3135388 and the MHC region are shown in Supplementary Tables 6 to 10. When excluding rs3135388, we observed an attenuation of the observed associations of the PRS; however, they remained statistically significant. We observed a larger decrease in the strength of the association after excluding the total MHC region, but the directions of the effect remained consistent (positive on global FA and negative on global MD).

To investigate tract‐specific associations between MS PRS and WM, the PRS with the strongest association with global FA (PT < 0.01) was regressed on individual WM tracts. Because the PRS showed no associations with global MD, we did not perform tract‐specific analyses for MD. The PRS showed positive associations with FA in 6 tracts: the superior longitudinal fasciculus (SLF left: β = 0.104, SE = 0.030, ΔR 2 = 0.0107, p = 5.90 × 10−4; SLF right: β = 0.080, SE = 0.030, ΔR 2= 0.0064, p = 8.17 × 10−3), the forceps minor (FMI; β = 0.091, SE = 0.030, ΔR 2= 0.0082, p = 2.75 × 10−3), the cingulate gyrus part of the cingulum (right: β = 0.069, SE = 0.029, ΔR 2= 0.0047, p = 0.02), and the corticospinal tract (CST left: β = 0.103, SE = 0.030, ΔR 2= 0.0104, p = 6.96 × 10−4; CST right: β = 0.072, SE = 0.030, ΔR 2= 0.0052, p = 0.02; Table 3, Fig 7). Exclusion of rs3135388 and the MHC region showed similar attenuations of the tract‐specific associations (see Fig 7).

Table 3

Regression Results of the Polygenic Risk Score (PT < 0.01) on the Fractional Anisotropy of 12 White‐Matter Tracts

Fractional Anisotropy, PT < 0.01
Tract	β	SE	ΔR 2	p
UNC L	0.048	0.030	0.0023	0.11
UNC R	0.048	0.030	0.0023	0.11
CGC L	0.032	0.029	0.0010	0.27
CGC R	0.069	0.029	0.0047	0.02b
SLF L	0.104	0.030	0.0107	5.90 × 10−4 a , d
SLF R	0.080	0.030	0.0064	8.17 × 10−3 a , c
FMI	0.091	0.030	0.0082	2.75 × 10−3 a , c
FMA	0.033	0.030	0.0011	0.28
ILF L	0.033	0.030	0.0011	0.27
ILF R	0.013	0.030	0.0002	0.67
CST L	0.103	0.030	0.0104	6.96 × 10−4 a , d
CST R	0.072	0.030	0.0052	0.02b

Regression results are corrected for age, sex, and 10 genetic principal components.

Significant p values after multiple testing correction (false discovery rate: p < 0.013).

Unrounded p < 0.05.

Unrounded p < 0.01.

Unrounded p < 0.001.

CGC = cingulate gyrus part of cingulum; CST = corticospinal tract; FMA = forceps major; FMI = forceps minor; ILF = inferior longitudinal fasciculus; L = left; R = right; SE = standard error; SLF = superior longitudinal fasciculus; UNC = uncinate fasciculus.

Figure 7

Effects of the polygenic risk score (p value threshold [PT] < 0.01 in our main study sample and PT < 0.005 in our replication sample) on the fractional anisotropy (FA) of 12 white‐matter tracts, including effects when removing rs3135388 and the major histocompatibility complex (MHC) region. Effects are standardized beta coefficients, corrected for age, sex, and 10 genetic principal components. An increased beta is represented by an increasingly yellow white matter tract; top to bottom: sagittal view, axial view, and coronal view. A = anterior; CGC = cingulate gyrus part of cingulum; CST = corticospinal tract; FMA = forceps major; FMI = forceps minor; ILF = inferior longitudinal fasciculus; L = left; P = posterior; R = right; SLF = superior longitudinal fasciculus; UNC = uncinate fasciculus.

In the significant tracts, the PRS of MS (PT < 0.01) was also regressed on RD and AD to explore the effect of the risk variants further (Supplementary Table 11). No significant associations were found for AD. For RD, the PRS had significant negative associations in the FMI (β = −0.078, SE = 0.030, ΔR 2= 0.0061, p = 9.92 × 10−3) and CST (left: β = −0.097, SE = 0.030, ΔR 2= 0.0094, p = 1.28 × 10−3; right: β = −0.081, SE = 0.030, ΔR 2= 0.0066, p = 7.07 × 10−3). When also correcting for global FA, none of the tract‐specific FA associations were significant (mean attenuation [%] = −0.37, range = −0.98 to 0.27).

Replication Sample

We tested our DTI findings for replication in an independent sample of 186 participants (median age = 8.5 years). Similar significant positive associations of the PRS of MS with global FA were observed (PT < 0.005: β = 0.189, SE = 0.072, ΔR 2= 0.0330, p = 9.40 × 10−3). We also found significant negative associations with global MD (PT < 0.005: β = −0.243, SE = 0.072, ΔR 2= 0.0546, p = 9.21 × 10−4; Supplementary Table 12).

Tract‐specific significant positive associations with FA were found in the forceps major, inferior longitudinal fasciculus, and the left CST (Supplementary Table 13, see Fig 7). Again, these associations attenuated to statistical nonsignificance after correction for global FA (mean attenuation [%] = −1.10, range = −3.40 to 1.05).

Inverse

To account for the relative overrepresentation of high maternal education in our DTI sample, we performed an additional analysis using IPCW based on maternal education. Although most of the effect estimates attenuated, the associations between the PRS for MS and global FA remained significant (Supplementary Tables 14–16).

Associations of

To investigate the effect of genetic risk for MS on cognitive functioning, we tested associations between the PRS and nonverbal IQ using a sample of 2,060 Caucasian participants (mean age at IQ assessment = 6.0 years, 50.3% female). The MS PRS was positively associated with nonverbal IQ across all thresholds (β = 0.057, SE = 0.022, ΔR 2= 0.0031, p = 0.01; Supplementary Table 17); however, these results did not remain significant following multiple testing correction. In addition, we performed an analysis using IPCW to account for possible differences in maternal education between the IQ study sample and the Generation R base sample. This showed the same consistent results (Supplementary Table 18).

In the participants who had DTI data available, in addition to the MS PRS and nonverbal IQ (n = 942; mean age at IQ assessment = 6.0 years, 49.0% female), the MS PRS was not significantly associated with nonverbal IQ across all thresholds (β = 0.052, SE = 0.033, ΔR 2= 0.0026, p = 0.12; Supplementary Table 19).

Discussion

We found evidence for a relationship between the combination of common genetic variants for MS and brain imaging phenotypes. Namely, we found a positive relationship between the MS PRS and WM microstructure in school‐age children from the general population. The MS PRS explained approximately 1% of the variance in global FA. Although the effects were small, we were able to replicate our findings in an independent sample of children in the Generation R Study who were scanned in an earlier wave and with a different MRI scanner. Effect directions in this sample were consistent although accompanied by increased beta coefficients compared to our main study sample, probably due to the smaller sample size having a greater chance of inflating the effect sizes.

Focusing on the specific WM tracts, we found significant positive associations with FA in the SLF, FMI, and CST. However, after correction for global FA, these relationships were no longer significant. This finding suggests a global effect of genetic MS risk on WM microstructure, instead of tract‐specific effects. We found no association between the genetic risk for MS and the global or regional T1‐weighted volumetric MRI measures.

In addition, we observed a positive, albeit nonsignificant, association between genetic risk for MS and nonverbal IQ, with the PRS explaining approximately 0.3% of its variance. This may be due to a higher global FA in children with higher genetic risk for MS.19 However, the temporal dynamics of our study (ie, global FA being measured relatively far after nonverbal IQ) stop us from carrying out a proper mediation analysis, which could provide insight into WM microstructure mediating the MS PRS and nonverbal IQ association.

Previous studies in MS patients report lower FA and higher diffusivity in several WM tracts of the CNS in adults and children.2, 4 Our findings in a population‐based study of children showed the opposite effect, providing an indication about the way the genetics of MS may influence WM and its components during neurodevelopment. Because FA and RD are associated with myelin integrity and myelination status of WM, our results may indicate a higher myelin integrity or different myelin status in children with a higher genetic risk for MS.38, 39 However, myelin is not the only factor determining FA. Branching and coherence of neuronal fibers also influence the degree of FA.40 A high number of crossing and branching fibers may reduce FA through an increased diffusion in multiple directions along the orientation of the fibers. With the degree of neuronal branching decreasing, the FA increases as the diffusivity along secondary and tertiary directions reduces (ie, RD).41 Hence, another explanation for the observed increased FA in children with a higher PRS accompanied by a lowered RD may be a significant decrease in the neuronal branching of several WM tracts. In this way, genetic risk of MS may affect WM maturation and disrupt neurodevelopment at the early age of our study sample. However, this mechanism would not be supported by the positive association we report between the MS PRS and nonverbal IQ. Future studies investigating the mechanisms involved in genetic MS risk, WM microstructure, and cognitive functioning would provide valuable insights. For example, an analysis exploring mediation by an altered WM microstructure on cognitive functioning would be appropriate.

Similar positive associations with the current study's diffusion metrics have been reported in previous studies investigating adults within the general population. A study by Ikram et al14 investigated the effect of MS risk variants in adults on WM microstructure and found similar directions of effect (a higher FA and lower diffusivity) for their top 5 most strongly associated risk variants when studying these diffusion parameters, although these findings did not survive multiple testing correction. A total of 110 non‐MHC risk variants were included, and the results may have been affected by age‐related brain changes. Brown et al15 found similar results (also not significant) in adults from the UK Biobank study when investigating associations between the PRS for MS and WM microstructure. Again, these findings could have been affected by age‐related brain changes and the broad age range of participants who were included in the study. Our pediatric cohort with a narrow age range minimizes the presence of age‐related differences; by using not only the genome‐wide significant variants in our genetic risk score, as done by Ikram et al,14 we attempted to reflect the polygenic architecture of MS as much as possible.

The association between our MS PRS, which predominantly consists of immunological variants, and structural brain findings seems counterintuitive. Susceptibility loci for MS are mainly expressed in different immune cells, including T cells and B cells, and microglia.9 Microglia in particular can affect oligodendrocyte functioning and myelination: through several interleukin pathways (ie, IL‐6 and IL‐1β), microglia are able to stimulate the survival and differentiation of oligodendrocyte precursor cells.42, 43, 44, 45 Because enrichment of MS risk variants has been shown in these pathways, oligodendrocytes and myelination could be influenced at a young age through alternative regulation of microglia function, leading to a differing WM microstructure.9

Considering the large contribution several MHC risk variants have in the genetic risk of MS, especially the class II HLA genes (for example HLA‐DPB1), we did not exclude MHC risk variants in the PRS of MS.46 Sensitivity analyses excluding these variants show a decrease in the magnitude of the effects, suggestive of a role of MHC risk variants in the different microstructure of WM. This may be driven by different MHC‐related gene expression in oligodendrocytes from individuals at a higher genetic risk for MS, which has been found in experimental autoimmune encephalomyelitis, an animal model mimicking the pathological process of MS.47

Taken together, individuals with a higher genetic risk for MS may show altered WM development via an interplay of altered immunological and oligodendrocyte function. The counterintuitive direction of effect (a positive relationship between MS PRS and FA, whereas other studies show a negative association in adult and pediatric MS patients) could indicate that environmental risk factors for MS play an additional role in disease pathophysiology and its WM alterations. How these environmental risk factors, such as smoking, vitamin D intake, or infection with the Epstein–Barr virus, may add to this process is unknown. Subsequently, how this altered WM integrity in individuals with a higher MS risk leads to an increased MS susceptibility remains elusive and yet to be investigated.

There are several strengths to this study. First, we studied a large sample of children drawn from the general population who underwent neuroimaging, giving us the possibility of investigating subtle differences between the MS PRS and the brain metrics. Second, all participants included in this study were scanned on the same study‐dedicated MRI scanner, eliminating possible interscanner differences from our main results. Third, our sample consisted of children from a narrow age range, which can reduce age‐related influences in WM development. Fourth, we replicated the findings in an independent group of children from the Generation R Study who were younger and scanned on a different MRI system.

Our study also has limitations: the association between the PRS for MS and imaging phenotypes was tested in a cross‐sectional study design, which precludes measuring changes over time. We did not include the first data collection wave of Generation R imaging participants in the main results of our study because these individuals were scanned on a different MRI system and encompassed a smaller subgroup.18 Instead, we used these participants as a replication sample to validate our results. Second, higher maternal education was overrepresented in our study. However, the height of the PRS was not affected by the level of maternal education, and after performing an IPCW analysis accounting for this overrepresentation of high maternal education, the association between genetic risk for MS and FA remained significant. Third, the GWAS for MS is relatively less powered compared to other GWAS studies, such as depression.48 This could lead to the low explained variance we report, indicating a subtle effect. Increasing the sample size of discovery GWAS studies can lead to an improved accuracy and power of studied associations.

Future studies are needed to explore the observed associations over time using a longitudinal study design. Additional insight into the genetic overlap between WM development and the pathophysiology of MS is also needed. For example, in addition to the effect of common variants, the association of rare MS risk variants with WM microstructure could be explored to provide a more complete understanding of disease pathogenesis.49 Furthermore, environmental risk factors of MS should be taken into account, as these are partially modifiable, to investigate the gene–environment interactions of MS and expand our knowledge about MS pathophysiology. In conclusion, we observed evidence that genetic risk for MS influences brain development, in particular the WM microstructure, in a pediatric population at an early age.

Author Contributions

C.L.d.M., P.R.J., T. J.W., and R.F.N. contributed to the conception and design of the study, and all authors contributed to the acquisition and analysis of the data and to drafting the text and preparing the figures.

Potential Conflicts of Interest

Nothing to report.

Table S1 Supplementary Material

Click here for additional data file.