In vivo marker of brainstem myelin is associated to quantitative sleep parameters in healthy young men

The regional integrity of brain subcortical structures has been implicated in sleep–wake regulation, however, their associations with sleep parameters remain largely unexplored. Here, we assessed association between quantitative Magnetic Resonance Imaging (qMRI)-derived marker of the myelin content of the brainstem and the variability in the sleep electrophysiology in a large sample of 18-to-31 years healthy young men (N = 321; ~ 22 years). Separate Generalized Additive Model for Location, Scale and Shape (GAMLSS) revealed that sleep onset latency and slow wave energy were significantly associated with MTsat estimates in the brainstem (pcorrected ≤ 0.03), with overall higher MTsat value associated with values reflecting better sleep quality. The association changed with age, however (MTsat-by-age interaction—pcorrected ≤ 0.03), with higher MTsat value linked to better values in the two sleep metrics in the younger individuals of our sample aged ~ 18 to 20 years. Similar associations were detected across different parts of the brainstem (pcorrected ≤ 0.03), suggesting that the overall maturation and integrity of the brainstem was associated with both sleep metrics. Our results suggest that myelination of the brainstem nuclei essential to regulation of sleep is associated with inter-individual differences in sleep characteristics during early adulthood. They may have implications for sleep disorders or neurological diseases related to myelin.

Sleep is essential to physical and mental health.Poor sleep is a predictor of mental health disorders, accelerated cognitive aging and neurodegeneration 1 .Although sleep undergoes profound changes over the lifetime 2 , it is remarkably stable within an individual over a shorter life period (e.g. a few months/years) 3 .The architecture of sleep is, however, highly variable across the individuals, even when they are very healthy 3 .A better understanding of such variability could provide important keys on the brain bases of a better sleep and for individually tailored sleep interventions.
The regulation of sleep heavily relies on nuclei of the reticular formation in the brainstem.They interact with each other and with nuclei of the diencephalon and basal forebrain to set vigilance state over the cortex 1,4 .Monoaminergic neurons, including the noradrenergic projections originating from the locus coeruleus (LC) and the serotonergic projection from the dorsal raphe nucleus, together with the cholinergic projections from laterodorsal tegmental nucleus (LDT) and the pedunculopontine tegmental (PPT) nucleus, and the dopaminergic projections of the ventral tegmental area (VTA) and substantia Nigra (SN) reach most cortical territories to promote wakefulness, when more active, or allow sleep when less active 4 .LC, LDT and PPT are also central to the switch between Rapid Eye Movement Sleep (REMS) and Slow Wave Sleep (SWS) as well as the oscillatory modes most typical of REMS and SWS [5][6][7] .The norepinephrine (NE)-LC system also contributes to the organisation of

Results
Following 3 weeks of stable sleep-wake timing, we recorded sleep of 321 healthy young men (22.1 years ± 2.7) under EEG (Table 1).To reduce the multiple comparison issue, we extracted only six sleep metrics; these were chosen because brainstem nuclei are reported to be involved in their regulation, at least based on animal research 4 , and because they cover our sleep characteristics of interest: (1) sleep onset latency (SOL), related to the initiation of sleep; (2) Slow wave energy (SWE) during SWS, corresponding to the cumulated overnight EEG power of delta band (0.5-4 Hz), an accepted marker of sleep need and SWS intensity 17 ; (3) sleep efficiency (SE; ratio between sleep time and time in bed), to assess overall sleep quality and continuity 18 ; (4) REMS percentage, to reflect the global architecture of sleep; (5) the cumulated theta power during REMS, associated with REMS intensity over its most typical oscillatory activity 19 and (6) number of arousals during REMS, to characterise sleep continuity 19 .SWE was preferred to slow wave activity, which consist in the EEG power of delta band over successive periods of time (typically a NREM-REM cycle for 1 h or 2 h), because it reflects the overall ability or need to generate slow wave and does not require modifications of the outputs of the validated automatic sleep stage scoring algorithm 20 through potentially subjective decisions (e.g. in case of skipped first REM episode).All participants subsequently underwent a qMRI protocol through which we computed the Magnetization Transfer saturation (MTsat) value, which is a semi-quantitative MRI measure linked to the myelin content 21 , for each voxel of the brain prior to extracting the average value of the bmGM compartment.The overview of the study design is provided in Fig. 1.
We first assessed whether our sleep metrics of interest varied with age.Despite the limited age range of our sample, we found significant decreases in SWE (p < 0.001), SE (p = 0.002), and REM theta power (p = 0.001) with age and significant increases in REM arousals (p = 0.001), which is in line with previously published agerelated changes in these variable [22][23][24] (Suppl.Table S1, Suppl.Fig. S1).No significant age-related changes in SOL (p = 0.11), and in REM percentage (p = 0.69) were detected.A negative borderline association in brainstem monoaminergic grey matter myelin content was further observed with age (p = 0.05) (Table S1, Suppl.Fig. S1).

Latency to sleep and SWS intensity are associated with the qMRI myelin marker over the monoaminergic brainstem compartment
Given the exploratory nature of the investigation, our statistical analyses consisted of GAMLSS which is a flexible distributional regression approach and is considered as an improvement and extension to the generalized linear models (GLM) and the generalized additive models (GAM) 25 .A first GAMLSS, with SOL as the dependent variable, yielded a significant negative main effect of the MTsat value of the bmGM (p = 2.2 × 10 −5 , p corr = 1.2 × 10 −4 ) and of age (p = 1.2 × 10 −4 , p corr = 7.2 × 10 −4 ) while controlling for body mass index (BMI), total sleep time (TST), and total intracranial volume (TIV), as well as MRI MPM sequence and scanner (see "Methods") (Table 2, Fig. 2A).The same GAMLSS also yielded an interaction between the bmGM MTsat value and age (p = 5.1 × 10 −5 , p corr = 0.0003).To gain insight in this interaction, we split our sample into 3 subsamples of similar size, respectively ranging from 18 to 20 years (N = 104), from 21 to 23 years (N = 133) and from 24 to 31 years    .Daytime sleepiness was measured by the Epworth Sleepiness Scale 61 , Chronotype was assessed by the Morningness-Eveningness Questionnaire (MEQ 88 ).Anxiety was estimated by the Beck Anxiety Inventory 58 and Mood was estimated by the 21-item Beck Depression Inventory II 59 ; Total Sleep Time (TST) was extracted from polysomnography recordings.a,b The data was available for 317 and 320 participants respectively.S2 for sleep metrics in each subgroup).We then computed GAMLSS for the 3 subgroups separately and found that while the association was significantly negative in the younger subsample (p = 0.006), i.e. with higher MTsat related to shorter SOL, there was no significant association with the intermediate and the older subsamples (p > 0.05) (Table 3, Fig. 2C).This suggests that the interaction between MTsat values and age that reflect a change in the slope and/or direction of the association between SOL and MTsat, is mostly driven by the younger individuals of our sample.In a second main GAMLSS, with SWE as the dependent variable, and controlling for the same factors as above, we found a significant positive main effect of the MTsat value of the bmGM (p = 0.01, p corr = 0.03) and of age (p = 0.04, p corr = 0.12) (Table 2, Fig. 2B).The GAMLSS also yielded an interaction between the bmGM MTsat value and age (p = 0.01, p corr = 0.03).As for SOL, to gain insight in this interaction, we split our sample into 3 subsamples.We did not find statistically significant association within a group (Table 3, Fig. 2D), and suspect that the interaction arises from the qualitative progressive switch in the direction of the association with a positive association in the younger subsample becoming negative in the older subsample.
Importantly, the main GAMLSS using the other sleep parameters of interest as dependent variables (sleep efficiency, REMS percentage, theta REMS power and No. of REMS arousals) were not significantly associated with the bmGM MTsat values (Table 2) suggesting that associations were specific to, or at least stronger for SOL and SWE.
Overall, these results show that increased MTsat values over the bmGM go along with enhanced sleep as reflected by a faster sleep onset and more intense NREM sleep.The association significantly changes from age 18-31 years, however.The younger individuals of our sample, aged 18-20 years, mirror the overall association with higher MTsat values associated with better sleep metrics while the association seem to progressively decrease or even revert in the intermediate subgroup, aged 21-23, and in the older subgroup, aged 24-31 years.

Associations between sleep metrics and the qMRI myelin marker is present across the different brainstem compartments
We assessed the regional specificity of the detected associations for the bmGM compartment and assessed potential associations between SOL and SWE with MTsat value computed over the other 2 brainstem compartments yielded by the segmentation procedure.The separate GAMLSS analyses with SOL or SWE, as dependent variable, and MTsat value in brainstem reticulate grey matter (brGM) and periaqueductal grey matter and posterior hypothalamus (bpGMpH) (while including age, BMI, TST, TIV, MPM sequence and scanner) yielded the same statistical outputs as the analyses focusing on the bmGM compartment.The associations between SOL and SWE with MT values are thus observed in other 2 brainstem compartments (Fig. 3A-D; Supplementary Fig. S2, Table S3).

Associations between sleep latency and intensity go beyond associations with the qMRI myelin marker over the prefrontal cortex
We further assessed the regional specificity of the finding for the brainstem by considering the MTsat value computed over the medial prefrontal cortex (mPFC), which shows progressive myelination during adolescence Table 3. Results derived from GAMLSS in age subgroups when testing for associations between sleep onset latency (SOL) or slow wave energy (SWE) and MTsat values computed over the brainstem monoaminergic gray matter compartment (bmGM).Significant values are in bold.BMI body mass index, bmGM brainstem monoaminergic gray matter compartment, MPM sequence type of multiparameter sequence, 2 sequence types were used (N = 222 for type 1; N = 99 for type 2-see "Methods"), MRI scanner data were acquired on 2 separate scanners (N = 286 for scanner 1; N = 35 for scanner 2-see "Methods"), MTsat magnetisation transfer saturation, SOL sleep onset latency, SWE slow wave energy, TIV total intracranial volume, TST total sleep time.www.nature.com/scientificreports/and early adulthood and is highly sensitive to sleep homeostasis 26,27 .Separate GAMLSS revealed significant main effects of the MTsat value over the mPFC for SOL (p = 0.01; p corr = 0.03) and SWE (p = 0.01; p corr = 0.03) (Suppl.Table S3; Fig. 3E,F) (while controlling for age, BMI, TST, TIV, MPM sequence and scanner).The GAMLSS using SOL and SWE as dependent variable yielded interactions between the mPFC MTsat value and age (respectively, p = 0.03 and p = 0.03; Table S3; Fig. S3).
Importantly, MTsat values over the mbGM and MTsat value over the mPFC were not correlated (Pearson's p = 0.17; Fig. S4).Critically, when including the mPFC MTsat values in the GAMLSS, the associations between SOL and SWE and the MTsat value computed over the bmGM remained significant both as main effect (SOL: p = 2.9 × 10 −5 , p corr = 1.7 × 10 −4 ; SWE: p = 0.01, p corr = 0.03) as well as in interaction effect between MT value and age (SOL: p = 6.7 × 10 −5 , p corr = 4.0 × 10 −4 ; SWE: p = 0.01, p corr = 0.03) (Suppl.Table S4).Therefore, when including the myelin marker over the prefrontal cortex, both SOL and SWE bear significant associations with MTsat value computed over the bmGM.As final validation of our main results, we computed cross-validation analyses by splitting dataset randomly into 70:30 train and test datasets which revealed that regressions between SOL or SWE and MTsat values computed over the bmGM remained significant in both trained and test dataset as compared to the when full dataset is used (Suppl.Table S5).

Discussion
The present study is part of the efforts to improve the understanding of the variability in healthy sleep.We tested whether qMRI, which informs about brain microstructural integrity 28 , would be correlated to canonical electrophysiological metrics of sleep, when estimated over the bmGM compartment which encompasses some of the most critical sleep-wake nuclei.We focussed on the MTsat metric because, in healthy tissue, it is a surrogate marker for myelin content 29,30 .We find that SOL and SWE are significantly related to brainstem MTsat values, with values reflecting better sleep composition linked to higher MTsat values.Despite the limited age-range of our sample, the associations changed with age, with MTsat values appearing mostly associated with better sleep composition in the younger subsamples.These associations are not specific to a particular compartment of the brainstem as they were also found when considering the other brainstem GM compartments.The associations with the brainstem are, however, explaining a different part of the variance in SOL and SWE than the significant association we also find for both metrics with MTsat values computed over the mPFC.Finally, none of the other  Myelin surrounds axons and has an intuitive impact on the crosstalk between distant brain regions reflected in the oscillations of the EEG 31 .Myelin is not only part of the white matter, but is also present in the grey matter, either at the bases of long axons or over short axons.Myelination can therefore also affect short-range neuronal connectivity and synchrony 31,32 .Myelination and synaptic pruning reflect a brain maturation process that progresses throughout childhood and adolescence and tail off through early adulthood, resulting in grey matter volume reduction and increase in white matter volume 14,15 .This maturation contributes in part to the important changes detected in sleep duration and quality taking place over the first 25 years of life 33 .There is for instance a significant reduction in the delta EEG frequency of SWS in adolescents that moves gradually with age from the posterior areas of the brain towards the anterior areas finishing at the PFC, a pattern that mirrors known brain myelination 34 and cannot be explained solely by a reduction in cortical thickness 35 .
Despite regional difference in myelin content across brainstem territories, an overall inverted U-shape reflect brainstem maturation with a progressive myelination up to ~ 25 years followed by a progressive decrease thereafter over the entire brainstem 16 .We find that several markers of sleep quality are associated with the brainstem myelin marker.While we expected specific correlations with the bmGM compartment, we find significant associations with all the three GM compartments of the brainstem.The fact that MTsat values over the brainstem and mPFC explain distinct parts of the variance in SOL and SWE, further support that regional maturation and integrity of the brainstem is important to some aspects of sleep.
The cross-sectional nature of our study does not allow inferring about causality.One can nevertheless posit that the association between MTsat and SOL or SWE are related to the maturation and integrity of the LC, raphe, SN, VTA, or cholinergic nucleus, which would optimize their local cross-talk-between neurons within and across nuclei-and their influence on cortical activity.This optimization would facilitate sleep onset and regulate slow wave production through and impact on the synchrony of neuron up and down states across patches of the brain.A relative silencing of the LC is for instance required to initiate sleep while transient noradrenalin release is locked to individual slow wave and spindles during SWS 4,36 .Similarly, activity of the serotoninergic neurons of the dorsal raphe (DR) is high during wakefulness, decreases during NREM sleep, and ceases during REM sleep 37 .Likewise, through their projection to the reticular nucleus of the thalamus, cholinergic neurons of the PPT and LDT facilitate sleep initiation and contribute to the production of spindle and therefore to the synchrony of cortical neuron activity 38,39 .The differences in myelination in the reticulate and PAG-posterior hypothalamus brainstem compartment could also contribute to the variability in the reported sleep indices.For instance, posterior hypothalamus on MRI images is likely to include lateral hypothalamic nuclei, which secretes orexin, a neuropeptide that regulates arousal, wakefulness, and appetite, and is involved in cataplexy 40 .
In vivo studies of adolescent animal models showed that disruption of slow wave sleep (N3 stage) results in significant alterations in the development of brain connectivity 14 .Still in rodents, oligodendrocyte precursor cells (OPC), programmed to generate myelin sheath and GM and WM 41,42 , proliferate two times faster during sleep than during wakefulness 43 .Sleep is also involved in synaptic pruning and long-term potentiation and depression of newly form neurites 44 , which may be wrapped in myelin sheets.The relationship between myelination and sleep has therefore been suggested to be bidirectional.In other words, our results are compatible with an impact of myelination of the diverse nuclei of the reticular formation on sleep quality as well as with an impact of sleep on myelination, where, for instance, repeated nights of sleep comprising more slow waves would increase myelin content.These hypotheses are not mutually exclusive and may take place concomitantly.Alternatively, our finding may also be incidental and arise from phenomenon that would both affect sleep and myelin.This phenomenon may also arguably be in part genetic which dictate brain myelin content and therefore affect sleep 26,45 .
Importantly, while higher myelin content, as indexed by MTsat, may favor or be favored by better sleep quality in early adulthood, we find that the association may be reversed after age ~ 21 years.Myelination follows a temporally symmetric time course across the adult life span and exhibit an inverted U-shape association with age, with a peak at ~ 25 years, in several brainstem substructures similar to what is observed in cerebrum, although substantially less pronounced 16 .We posit that our findings could mean that delay in the maturation of the brainstem across the ascending part of this inverted U-shape, i.e. over the second and third decade of life, is associated with lower sleep quality-with significant difference starting at ~ 20-25 years.Our interpretation is, however, based on significant interaction that we could partially draw from the significant association with specific age subgroups.Our findings are in line with the previously proposed occurrence of bi-directional interactions between sleep and myelin plasticity and endorse that sleep-wake-plasticity interactions may occur not only in age-specific manner but also in brain region-specific manner 46 .In line with this assumption, we recently reported in a whole-brain analysis of large part of the current sample that early night frontal SWA was negatively associated with regionally decreased myelin estimates, as indexed through MTsat values, in the temporal portion of the inferior longitudinal fasciculus 47 .These results and the present ones show that regional myelin change may affect different aspects of sleep.A full appreciation of these region-specific associations over the lifespan requires more investigation, for instance using a sample with a larger age range of the healthy participants.
Our study bears limitations.We only included men and did not include a longitudinal component to the protocol.Calibration of the actimetry devices was not performed prior to each recording so proper measures of physical activity could not be derived and included in any of the analyses.To circumscribe the multiple comparison issue, we also only focussed on six sleep metrics and a single qMRI parameter.In addition, MTsat values are not specific to myelin but rather correlates to all macromolecules of the tissue.In healthy individuals, this macromolecular content mainly concerns cell membrane and therefore myelin of oligodendrocytes 48 .MTsat may also be more directly related to the density of neurons and therefore indirectly related to myelin 49 .Other qMRI parameters depend to a lesser extent on myelin such as R2*, which is mainly driven by iron content 50 .
Vol:.( 1234567890 It would be of interest to consider other qMRI quantification approaches along with MTsat as it would be to focus on other parts of the brain.Scarce studies reported that changes in the volumes of the hypothalamus and thalamus were associated to alterations in sleep, though mostly in clinical populations (narcolepsy, obstructive sleep apnea, or neurodegeneration) 51,52 .Also, a reduction in thalamic GM density which could arguably result from synapses and myelin loss was reported to mediate part of the oscillatory changes in the EEG with ageing 13 .Finally, several nuclei of the brainstem are keys to sleep initiation and regulation and could contribute to the associations we detect 4 .Other type of data (higher resolution, other contrasts, etc.) and/or the development of other segmentation tools (e.g. 53) are required to isolate specific nuclei.
Our study contributes to a better understanding of the brain correlates underlying sleep variability in healthy individuals.It may also have implication for patient populations that suffer from sleep disorders or degradation of myelin or both.For example, lesions at the brainstem and spinal cord level are more frequently associated with the appearance of REMS behaviour disorder (RBD) and restless legs syndrome (RLS) which constitute risk factors for Parkinson's disease 54 .Likewise, sleep complaints are frequent in multiple sclerosis 55,56 and, based on our finding, myelin degradation over the brainstem could contribute to such complaints.

Ethics statement
Approval for this study was obtained from the Ethics Committee of the Faculty of Medicine at the University of Liège, Belgium.Written informed consent was obtained from each participant prior to their participation and they were financially compensated.The study was performed in accordance with relevant institutional guidelines/ regulations and also in accordance with the Declaration of Helsinki.

Participants
The study sample and study protocol is as previously published 57 .We recruited a cohort comprising 364 healthy young men (aged 18-31) recruited between 2011 and 2015 as part of a large study assessing association between sleep architecture and circadian rhythmicity incorporating genetics and polygenic risk score assessments (for details see 57 ).The study only incorporated men to increase genetic homogeneity.The participants were excluded based on the following criteria: body mass index (BMI) > 27; presence of psychiatric history or severe brain injury; addiction; chronic medication affecting the central nervous system (CNS); smoking, excessive alcohol (> 14 units/week) or caffeine (> 3 cups/day) intake; shift work in the past year; trans-meridian travels in the past 3 months; presence of moderate to severe subjective anxiety and depression as measured by the Beck Anxiety Inventory (BAI; score > 16) 58 and Beck Depression Inventory II (BDI; score > 19) 59 , respectively.Participants that exhibited poor sleep quality as assessed by the Pittsburgh Sleep Quality Index (PSQI) 60 (score > 7), excessive daytime sleepiness as index by the Epworth Sleepiness Scale 61 (ESS; score > 14) or excessive sleep apnea (apnea-hypopnea index > 15/h; 2017 American Academy of Sleep Medicine criteria, version 2.4) were excluded based on an in-lab screening night of polysomnography.
Out of 364, 28 participants were excluded after quality assessment of MRI data (12 had incomplete MRI acquisition, 15 had hyperintensity issues or movement artifacts in the MRI images, 1 had problem in segmentation process), 14 participants had incomplete baseline sleep electrophysiological data and 1 outlier resulting in a final sample of 321 participants.The characteristics of the final participant cohort are reported in Table 1.

Sleep protocol
As described in [Muto 2021], individual sleep-wake history was strictly controlled: during the 3 weeks preceding the in-lab experiment, participants were instructed to follow a regular sleep schedule according to their habitual sleep timing (± 30 min for the first 2 weeks; ± 15 min for the last week; verified using actigraphy data-Actiwatch 4, CamNtech, Cambridge, UK).A urine drug test was performed (Multipanel Drug test, SureScreen Diagnostics Ltd) before completing an adaptation night at habitual sleep/wake schedule during which a full polysomnography was recorded to screen for sleep-related breathing disorders or periodic limb movements.On day 2, participants left the lab in the morning with the instruction not to nap which was confirmed with actigraphy data.They came back to the laboratory at the end of day 2 (3.5 h hours before scheduled lights-off) and completed a baseline night of sleep (BAS) in full darkness centered on the average sleep midpoint of the preceding week.The current study focuses only on the baseline night of sleep.The remaining of the protocol included successively an extended nighttime sleep opportunity period, a daytime nap period, a normal night, a 40 h sleep deprivation under constant routine condition (light < 5 lx) and finally a 12 h recovery night.Outside sleep opportunity and sleep deprivation periods, participants were maintained in normal room light levels ranging between 50 and 1000 lx depending on location and gaze.The variability in the actual light received due to gaze position and movement was not documented.

EEG acquisitions and analyses
Polysomnographic sleep data were acquired using Vamp amplifiers (Brain Products, Germany).The electrode montage consisted of 10 EEG channels (F3, Fz, F4, C3, Cz, C4, Pz, O1, O2, and A1; reference to right mastoid), 2 bipolar EOGs, 2 bipolar EMGs, and 2 bipolar ECGs.Acquisition on the screening night of sleep also included respiration belts, oximeter and nasal flow, 2 electrodes on one leg, but included only Fz, C3, Cz, Pz, Oz, and A1 channels.EEG data were re-referenced off-line to average mastoids.Scoring of sleep stages was performed using a validated automatic algorithm (ASEEGA, PHYSIP, Paris, France) in 30-s epochs 62 and according to 2017 American Academy of Sleep Medicine criteria, version 2.4.An automatic artefact and arousal detection algorithm with adapting thresholds 63 was further applied and artefact and arousal periods were excluded from subsequent analyses.Power spectrum was computed for each channel using a Fourier transform on successive 4-s bins,  64 .Thus we computed slow wave energy (SWE)-cumulated power in the delta frequency band during N2 and N3 sleep stages, an accepted measure of sleep need 65 , and similar to that we computed the cumulated theta (4-8 Hz) power in REM sleep.We then computed the cumulated power over the remaining EEG bands, separately for NREM and REM sleep: alpha (8-12 Hz), sigma (12-16 Hz), beta (16-25 Hz) and theta (4-8 Hz) bands.As the frontal regions are most sensitive to sleep-wake history 66 , we considered only the frontal electrodes (mean over F3, Fz, and F4), as well as to facilitate interpretation of future large-scale studies using headband EEG, often restricted to frontal electrodes.Our analyses focused on six sleep metrics to limit issues of multiple comparisons while spanning the most important aspects of sleep EEG:

Image acquisition and quantitative maps creation
On average 3 weeks upon completion of the sleep study protocol, participants were examined on a 3 T MR system.Out of the 321, neuroimaging data of 286 participants were acquired on a 3 T head-only MRI-scanner (Magnetom Allegra, Siemens, Erlangen, Germany).Further, out of these 286 participants, 187 were scanned using an Echo planar imaging (EPI) sequence, while 99 were scanned using an Actual Flip Angle Imaging (AFI) sequence.MRI data for the remaining 35 participants were acquired on a 3 T whole-body MRI-scanner (Magnetom Prisma, Siemens, Erlangen, Germany), with EPI sequence due to scanner replacement.The acquisition included a whole brain quantitative multiparameter protocol (MPM) as described in 67 and 68 .
The MPM protocol used in the present study has been already validated for multi-centric acquisitions 28 .Briefly, it consists of three co-localized series of 3D multi-echo fast low angle shot (FLASH) acquisitions at 1 × 1 × 1 mm 3 resolution and two additional calibration sequences to correct for inhomogeneities in the Radio Frequency (RF) transmit field 69 .The FLASH data sets were acquired with predominantly proton density (PD), T1, and magnetisation transfer (MT) weighting, referred to in the following as PDw, T1w and MTw echoes.Volumes were acquired in 176 sagittal slices using a 256 × 224 voxel matrix.
Quantitative multi-parametric volumes (PDw, T1w and MTw) were auto-reoriented against an MNI template and MT saturation (MTsat), PD, R1 and R2* maps were created using the hMRI toolbox 67 (http:// hmri.info) implemented as an add-on toolbox for SPM12 (Statistical Para-metric Mapping, Wellcome Centre for Human Neuroimaging, London, UK, http:// www.fil.ion.ucl.ac.uk/ spm, University College London, revision 12.4).Here, we focused on the magnetization transfer saturation map (MT) which is related to the exchange of magnetization between mobile water protons that are bound to macromolecules as found in myelin.MT intensity values have been closely related to myelin content as shown in postmortem studies 21 .Contrary to the commonly used MT ratio (percentage reduction in steady state signal), the MT saturation map explicitly accounts for spatially varying T1 relaxation time and flip angles 70 .The map creation module includes the determination of B1 transmit bias field maps for transmit bias correction of the quantitative data.For this dataset, two different methods were used, the EPI 69 and Actual Flip Angle Imaging (AFI) method 71 , according to the sequence implemented at the time of acquisition.

Quantitative multiparametric MRI data analysis
For the extraction of quantitative MTsat-derived myelin markers from specific brainstem regions MT and PD maps were first segmented into grey, white, and CSF tissue class maps using Unified Segmentation (US) within SPM12 72 .The whole-brain segmentation outputs were then diffeomorphically registered to a study-specific template, compatible with the MNI space, created using Shoot toolbox in SPM12 73 in order to generate deformation fields that were used to warp MT and PD maps into the study-specific average space.Brainstem segmentation was then performed using unified segmentation with brainstem subregion tissue probability maps that were generated according to previously described methods based on a modified multivariate mixture of Gaussians 74 .Out of the 4 brainstem tissue classes produced, only tissue class 1 was considered for our main analysis, as it includes the substantia nigra, locus coeruleus and raphe nuclei designated as "brainstem monoaminergic grey matter, " (bmGM).The other two tissue classes mainly contained dorsal cranial nerve nuclei and were used for exploratory analyses only [tissue class 2: nucleus reticularis throughout its length and the pontine nuclei (reticulated grey matter; brGM); tissue class 3: Periaqueductal grey matter and posterior hypothalamus (bpGMpH) 75 .The tissue class 4: brainstem white matter (bWM) was excluded from the analysis as white matter tracts cross the brainstem and LC projections are mainly unmyelinated.The 3 tissue classes were warped back to individual space, using inverse deformation fields.Tissue class-specific smoothing (full width at half-maximum of 3 mm isotropic) was applied, and mean tissue class images were created after averaging in SPM Masking toolbox.
For MTsat-derived myelin markers from mPFC, segmented images were normalized using geodesic shooting 76 and smoothed with a tissue-weighted, for GM and WM separately, kernel of 4 mm FWHM.Using the TD-ICBM Human atlas as implemented in the WFU-Pickatlas toolbox v 3.0.5bfor SPM12 77 we defined a mask the medial PFC via TD Brodmann areas in MNI space.For creating medial PFC mask, BA 24, BA 25 and BA 32 regions were merged.The modulated spatially normalised tissue maps of grey matter warped into MNI space were used for voxel-based quantification analyses (VBQ).Matlab based REX toolbox (https:// web.mit.edu/ swg/ softw are.htm) was used to extract the mean MTsat intensity values (single-subject beta values) across each ROI mask.

Statistical analysis
All analysis was carried out within R environment (version > 4.1.3)(R Development Core Team, 2017).We employed generalized additive models for location scale and shape (GAMLSS) 78,79 to separately test the associations between six sleep metrics of interests (SOL, SWE, Sleep efficiency, REM percentage, REM power and No. of arousals in REM), as dependent variable, and the estimated MTsat values from bmGM as an independent variable.GAMLSS was selected based on data distribution and lowest AIC values among General Linear Models (GLM) and GAM models.In order to control for potential association between the sleep metrics and cortical MTsat values, average MTsat value was computed over the medial prefrontal cortex, as it is located below the electrodes of interest and an important source of slow waves 80 .We further checked the specificity of associations with bmGM tissue class and computed MTsat intensity values from other brainstem tissue classes (brGM and bpGMpH) used as independent variable for subsequent analysis.bWM was not included in the analysis as it may be extremely heterogeneous, with motor (pyramidal) and somatosensory tracts.GAMLSS are univariate distributional regression models, where all the parameters of the assumed distribution for the dependent variable can be modelled as additive functions of the predictor/independent variables.On the other hand, GLMM and GAM are restricted to the exponential family of distributions 78 .GAMLSS offer a large variety of distributions with up to four parameters-location (e.g.mean), scale (e.g.variance), shape (e.g.skewness) and shape (e.g.kurtosis), classically noted as µ, σ, ν and τ.While only µ is modelled in (G)LM(M) and GAM(M), in GAMLSS all four parameters can be modelled, either with linear parametric, non-linear parametric or non-parametric (smooth) functions of the predictors.GAMLSS algorithms have been mainly designed to complete two tasks: maximize a penalized log-likelihood function addressing the estimates of fixed and random parameters, and evaluate the various smoothing parameters appropriately [81][82][83] .
In the present study, either generalised Gamma (GA) or Gamma Lopatatsidis-Green (GG) family distribution was used for the GAMLSS regression modelling based on the Q-Q plot.Age, BMI, total sleep time (TST) and total intracranial volume (TIV) were included as covariates.Further, MPM sequence type and MRI scanner type was also included in the regression models to control for the protocol difference and scanner type during acquisitions.For age-subgroup interaction analysis plots GAM was used with either SWE or SOL as DV and age as categorical independent variable (IV) x MT value as continuous IV.Prior to the analysis, influential outliers were screened using worm plot.The worm plot (a detrended Q-Q plot) is a diagnostic tool for checking model fit, improving the fit and comparing the fit of different models 84 .Cross-validation analyses were used for validating the models by splitting dataset randomly into 70:30 train and test datasets and regression performance was measured using the Root Mean Square Error (RMSE) metric.Because of the exploratory nature of the hypotheses, Benjamini and Hochberg False Discovery Rate (FDR) correction for 6 independent tests was used to test for significant associations.
Optimal sensitivity and power analyses in GAMLSS remain under investigation 85 .We nevertheless computed a prior sensitivity analysis to get an indication of the minimum detectable effect size in our main analyses given our sample size.According to G*Power 3 (version 3.1.9.4) 86,87 taking into account a power of 0.8, an error rate α of 0.01 (corrected for 5 tests), a sample size of 321 allowed us to detect small effect sizes r > 0.24 (2-sided; absolute values; confidence interval: 0.13-0.34;R 2 > 0.06, R 2 confidence interval: 0.017-0.12)within a linear multiple regression framework including 2 tested predictor (MT value, age) and 5 other covariates (BMI, TST, TIV, MPM sequence, Scanner type).

Figure 1 .
Figure 1.Overview of the study design: (A) In-lab recordings of habitual sleep to extract sleep macro and microstructure metrics.(B) 7 T MRI multiparameter protocol (MPM)acquisitions to generate quantitative (qMRI) maps.Magnetization transfer (MT) value (related to myelin content) was averaged was computed over the brainstem monoaminergic gray matter compartment (bmGM) compartment.PDw proton density weighted, MTw magnetization transfer weighted, T1w T1 weighted.

Figure 2 .
Figure 2. Plots (A,B) depict Pearson correlations between sleep onset latency (SOL) and slow wave energy (SWE) with MTsat intensity in bmGM respectively (N = 321), see Table 2 for statistical outputs of GAMLSSs.Plots (C,D) shows association between sleep parameters and MTsat values by age group for SOL and SWE with MTsat intensity in brainstem monoaminergic grey matter (bmGM) respectively, see Tables2 and 3for statistical outputs of GAMLSSs.

.
60aracteristics of the final participant cohort included in the analyses.Sleep quality was assessed by the Pittsburgh Sleep Quality index (PSQI60

Table 2 .
Results derived from GAMLSS when testing for associations between sleep parameters and MTsat values computed over brainstem monoaminergic grey matter (bmGM) with age as interacting variable.Significant values are in bold.BMI body mass index, bmGM brainstem monoaminergic gray matter compartment, MPM sequence type of multiparameter sequence, 2 sequence types were used (N = 222 for type 1; N = 99 for type 2-see "Methods"), MRI scanner data were acquired on 2 separate scanners (N = 286 for scanner 1; N = 35 for scanner 2-see "Methods"), MTsat magnetisation transfer saturation, REMS rapid eye movement sleep, SE sleep efficiency, SOL sleep onset latency, SWE slow wave energy, TIV total intracranial volume, TST total sleep time.

Table 2
for statistical outputs of GAMLSSs.Plots (C,D) shows association between sleep parameters and MTsat values by age group for SOL and SWE with MTsat intensity in brainstem monoaminergic grey matter (bmGM) respectively, see Tables2 and 3for statistical outputs of GAMLSSs. ) Vol:.(1234567890) Scientific Reports | (2023) 13:20873 | https://doi.org/10.1038/s41598-023-47753-x TableS4for statistical outputs of GAMLSSs.metrics of interest, related to overall SWS-REMS organisation and to REMS fragmentation and intensity, were found to be significantly correlated to MTsat value over the brainstem, suggesting that associations are stronger for, or specific to SOL and SWE.