Prader-Willi syndrome (PWS) is a rare genetic neurodevelopmental disorder caused by the loss of paternal expression of maternally imprinted genes on chromosome 15q11-13 [1]. The genetic subtypes of PWS include paternal deletion, maternal uniparental disomy and imprinting center defects [2]. The trajectory of the PWS phenotype primarily ranges from severe hypotonia in infancy to an insatiable appetite during the preschool period, with the appearance of hormone deficiencies, behavioral problems and dysautonomia [3]. The clinical features of PWS are characterized by short stature, bulimia, irritability, skin scratching, dysmorphic facial features, high threshold and developmental delays [4]. These developmental, mental and behavioral issues inevitably lead to structural and functional changes in the brains of PWS patients [5,6,7,8].

Current multimodal magnetic resonance imaging (MRI) has been used to assess brain alterations in structural and functional domains in PWS [9]. Rice et al. [6] showed white matter lesions in the internal capsule and the splenium of the corpus callosum in young adults with PWS; these findings suggested that white matter microstructural abnormalities were linked to cognitive and emotion-related functions, such as task switching, emotion recognition, semantic processing, and sensorimotor deficits. Zhang et al. [7] selected four neural networks and discovered altered functional connectivity among these brain regions; their findings provided additional evidence for an obesity model of overeating and even addiction. Our recent study focused on explaining specific developmental delays in children with PWS at the level of large-scale functional networks based on resting-state functional MRI (rs-fMRI) [8]. Until now, none of these studies have directly considered combining structural and functional networks in PWS during early childhood.

Structural-functional coupling has recently been proposed to combine neural structural and functional networks through quantitative MRI studies [10]. This coupling pattern may detect subtle shifts in the associations between structural networks (physical metabolic substrates) and functional networks (physiological metabolic activity) in disease. For instance, Chen et al. [11] reported increased structural-functional coupling strength in patients with idiopathic tinnitus for whom sound therapy was ineffective. Altered structural-functional coupling implies inconsistencies at the macroscopic level of structural and functional networks, and further quantitative analysis of complex brain networks is necessary to reveal subtle brain pathophysiological abnormalities. Topological network analysis based on graph theory can further evaluate global and nodal parallel information processing in complex brain networks. Indeed, in common neuropsychiatric disorders, specific connectivity patterns may be responsible for the global and nodal topological properties of brain networks [12,13,14]. However, network-based studies examining alterations in structural-functional coupling (at the macro level) and intrinsic topological properties (at the micro level from global to nodal) in children with PWS are rare.

This study aimed to investigate large-scale network alterations associated with early developmental delays in children with PWS in a stepwise manner using diffusion tensor imaging (DTI) and rs-fMRI data. First, we constructed brain structural and functional networks and hypothesized that decoupling of these two modalities exists at the macroscopic level in children with PWS. Second, we further focused on the distinct changes in the global and nodal topological properties of these two modalities in children with PWS. Finally, we examined the relationships between brain network-based changes and PWS-related developmental states.

Materials and methodsParticipants

The Human Research Ethics Committee of Children’s Hospital of Chongqing Medical University approved this protocol. Written informed consent was obtained from the parents of all subjects in accordance with the ethical guidelines of the Declaration of Helsinki [15]. The registration number on ClinicalTrials.gov is ChiCTR2100046551.

Thirty children with PWS and 32 healthy controls (HCs) were enrolled in this prospective observational study. Children with PWS were recruited from the Endocrinology Department of Children’s Hospital of Chongqing Medical University. PWS was diagnosed according to the guideline criteria in 2011 [4]. All HCs with typical development were recruited from a population undergoing clinical MRI scans of non-neurological body regions, ensuring that all participants demonstrated an absence of any notable neurological developmental deviations. None of the individuals had a history of neurological or psychiatric disorders, nor were they taking any kind of medication. Five patients and 4 HCs were excluded due to large head movements (> 2 mm translation or > 2° rotation) during scanning. The final data set contained information on 25 PWS patients (42.52 ± 22.55 months) and 28 HCs (44.55 ± 16.37 months) (Table 1). All the subjects underwent laboratory tests and MRI scanning. Subjects with PWS all underwent the Griffiths Development Scale (GDS) assessment before MRI scanning. Developmental delays were noted when the general quotient or subscale quotient was at least 2 standard deviations below the mean [16].

Table 1 Demographic, clinical and developmental assessment results of participantsMRI data acquisition

MRI data were collected with a 3.0 T Philips Achieva scanner with an 8-channel head coil. Foam padding was used to reduce head motion, and earplugs were used to reduce scanner noise. To improve the quality of neuroimages, all participants were under moderate sleep deprivation achieved with intravenous administration of propofol (loading dose of 1 mg/kg, followed by 200–300 µg/kg/min).

DTI data were acquired using spin echo-based echo planar imaging sequences with the following parameters: repetition time (TR) = 7865 ms, echo time (TE) = 66 ms, slice thickness = 2 mm, gap = 72 mm, slices = 70, field of view (FOV) = 224 mm × 224 mm, acquisition matrix = 112 × 112, flip angle = 90°, voxel size = 2 mm × 2 mm × 2 mm, nonzero b value = 1000 s/mm2, and total time = 5 min 41 s.

High-resolution T1-weighted imaging data were acquired using a three-dimensional turbo fast echo sequence, and the parameters were as follows: TR = 7.4 ms, TE = 3.8 ms, slice thickness = 1 mm, gap = 0 mm, slices = 260, FOV = 250 mm × 250 mm, acquisition matrix = 228 × 227, flip angle = 8°, voxel size = 0.60 mm × 1.04 mm × 1.04 mm, and total time = 4 min 16 s.

Rs-fMRI data were acquired using a gradient echo-planar imaging sequence, and the parameters were as follows: TR = 2000 ms, TE = 35 ms, slice thickness = 4 mm, slices = 33, FOV = 240 mm × 240 mm, acquisition matrix = 80 × 78, flip angle = 90°, voxel size = 3.75 mm × 3.75 mm × 4 mm with no gap, and total time = 8 min 06 s.

Preprocessing and network construction

Structural preprocessing steps (DTI data) were conducted by using the pipeline toolbox PANDA [17] based on FSL software (https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/). The resolution of each voxel was resampled to 2 mm × 2 mm × 2 mm, and the raw images were cropped. The B0 images (b = 0 s/mm2) were acquired for coregistration using affine transformations, with the eddy current distortions and head motions corrected. Then, whole-brain fiber tracking was performed using the Fibre Assignment by Continuous Tracking algorithm. We chose a fractional anisotropy (FA) of 0.15 or a turning angle greater than 45° as the continuous fiber termination criterion [11].

Functional preprocessing steps (rs-fMRI data) were conducted by using DPABI (http://www.rfmri.org/dpabi) [18] in MATLAB 2013b. The first 10 volumes were discarded from each functional time series. Slice-timing correction and realignment for head motion correction were performed. Subjects with head motion greater than 2 mm in translation or 2° in rotation were excluded. The corrected data set was normalized to the Montreal Neurological Institute (MNI) template with a resolution of 3 mm × 3 mm × 3 mm. Then, the covariates were regressed from the time course of each voxel, including head motion parameters in 24 directions and four average confounding signals (cerebrospinal fluid, white matter, gray matter and the whole brain). Finally, the time series were temporally bandpass filtered (0.01–0.08 Hz).

Structural and functional networks were obtained through the GRETNA toolbox (http://www.nitrc.org/projects/gretna/) [19] by using the automated anatomical labeling (AAL) algorithm. The nodes of the structural/functional networks were defined based on 90 noncerebellar anatomical regions of interest from the MNI space [20]. For the structural networks, the edges were obtained by the averaged FA of the linking fibers between every pair of node weights. For functional networks, the edges were acquired by computing Pearson correlation coefficients between the processed temporal series of every pair of nodes. For the structural/functional network construction flowchart details, see Fig. 1.

Structural-functional coupling analysis

Pearson’s correlation was calculated to quantify the strength of structural-functional coupling for each subject [20]. Before that, the nonzero structural network edges were extracted, rescaled to a Gaussian distribution, and then related to the corresponding functional network edges.

Global and nodal network topological property calculations

To evaluate the topological properties of structural and functional networks, both global and nodal graph theoretical metrics were calculated by using the GRETNA toolbox. We used a matrix element value (threshold sequence set to 0) for structural networks and a connection sparsity threshold (from 0.05 to 0.5, with an interval of 0.01) for functional networks. We further calculated the area under the curve (AUC) over the sparsity for intergroup differences, providing a summarized scale for the topological characterization of brain networks independent of a single threshold selection. For global properties, the included metrics were the normalized clustering coefficient (γ), normalized characteristic path (λ), small-worldness (σ; = γ/λ), clustering coefficient (Cp), characteristic path length (Lp), global efficiency (Eg) and local efficiency. For nodal properties, betweenness centrality (Bc), degree centrality (Dc) and nodal efficiency (Ne) were applied. For the detailed calculation methods, please see the Supplementary material.

Statistical analysis

Statistical analyses were performed by using Statistical Product and Service Solutions version 23.0 (SPSS, IBM Corp., NY, United States). The Kolmogorov‒Smirnov test was performed on demographic and laboratory data to test the normality of continuous variables. Intergroup differences in normally distributed and non-normally distributed continuous data were assessed by using two-tailed independent samples t tests and Mann‒Whitney U tests, respectively. Intergroup differences in the categorical data (such as sex) were examined by using the χ-squared test. The results are reported at the significance level of p < 0.05.

Intergroup differences in structural-functional coupling, global and nodal properties were analyzed by using general linear models (GLMs) on graph theoretical measures controlling for age and sex. Structural – functional coupling and global properties were tested by using permutation tests. This procedure was repeated for 5000 permutations, setting a threshold at p < 0.05. For global properties, p < 0.05 was set as the significance level. For nodal properties, multiple comparisons were corrected using a false-positive correction p < 1/90 = 0.011 [21].

Partial correlation analyses were used to identify possible relationships between network features (structural-functional coupling and global/nodal properties) and GDS in children with PWS with age and sex as nuisance covariates, setting a threshold at p < 0.05.