Introduction

The COVID-19 pandemic has presented unprecedented challenges to both the physical and psychological well-being of individuals worldwide1. Now, in the fifth year of this global health crisis, mounting evidence indicates that SARS-CoV-2 exerts a detrimental impact on the central nervous system, particularly the brain. This viral influence has been linked to a range of neurological conditions, including headaches, dizziness, anosmia, ageusia, and cognitive impairments, all of which point to the virus’s potential neurotropic properties2,3,4. Such findings underscore the virus’s capacity to infiltrate the nervous system, posing significant challenges in understanding and managing its long-term neurological consequences5.

Magnetic Resonance Imaging (MRI) has become a critical non-invasive technique for investigating the neuropathological mechanisms underlying brain alterations associated with COVID-196,7. Several MRI studies have reported significant reductions in brain size, thinning of grey matter, disruptions in functional connectivity, and bilateral declines in glymphatic function in COVID-19 patients—changes that remain observable even months after recovery7,8,9,10,11,12. These structural and functional alterations may underlie the persistent neurological symptoms that many COVID-19 survivors continue to experience13.

The early years of childhood are characterized by rapid, dynamic changes in brain structure, behavior, and cognition14,15. At the microstructural level, neurogenesis and synaptogenesis are occurring at a rapid pace, while the blood-brain barrier remains underdeveloped. At the macroscopic level, the brain undergoes significant morphological changes that further refine the separation and integration of brain functions16,17,18. Given that the first few years of life represent a critical window for brain development, it is plausible to hypothesize that COVID-19 infection could similarly affect brain structure and function in young children. While numerous studies have explored the impact of COVID-19 on the adult and elderly brain, research focusing on pediatric populations remains limited, primarily consisting of case reports and meta-analyses19,20. Consequently, there is a lack of quantitative data regarding the virus’s impact on the developing brain in children.

To address this critical gap, we recruited a cohort of children who had contracted mild COVID-19 and compared them with healthy controls (HCs). High-resolution structural and diffusion MRI data were collected to assess alterations in brain structure and function. Whole-brain vertex-wise morphometric analyses were performed, and structural covariance networks (SCNs) were constructed to evaluate the effects of COVID-19 on brain structure. SCNs are valuable for understanding how brain development, disease states, and shared environmental factors influence brain morphology21,22. Furthermore, diffusion tensor imaging along the perivascular space (DTI-ALPS), which assesses fluid flow along perivascular spaces crucial for waste clearance, and choroid plexus (CP) volume, which evaluates potential alterations in cerebrospinal fluid (CSF) production, were employed as MRI proxies to examine the function of the glymphatic system23.

Methods

Participants

This study was conducted in accordance with the Declaration of Helsinki (revised in 2013) and received approval from the Ethics Committee of the Children’s Hospital of Fudan University (Xiamen Branch), Xiamen Children’s Hospital. Informed consent was obtained from the parents or guardians of all participants after a thorough explanation of the study’s objectives. A total of 19 children (16 males, 3 females; mean age 3.4 ± 1.9 years) diagnosed with mild COVID-19 were enrolled based on the following inclusion criteria: (1) confirmed SARS-CoV-2 infection via antigen or polymerase chain reaction (PCR) testing; (2) at least two weeks post-recovery from infection, defined as both a negative SARS-CoV-2 test result and resolution of symptoms; (3) aged 1–6 years; (4) no history of primary or secondary neurological disease, such as craniocerebral injury, tumors, epilepsy, or meningitis; (5) no major or genetic conditions affecting growth, development, or cognition; and (6) no contraindications to MRI. The exclusion criteria were: (1) incomplete image acquisition; (2) poor image quality; and (3) inadequate image segmentation during preprocessing. The median interval between the onset of COVID-19 and MRI examination was 30.0 days (interquartile range [IQR]: 27.5, 34.5 days). Additionally, 22 age-comparable HCs (13 males, 9 females; mean age 3.9 ± 1.3 years) were recruited between February 2023 and May 2023. These controls had no history of SARS-CoV-2 infection (confirmed by antigen or PCR testing and self-reported history), no neurological disease, and no contraindications to MRI. All participants were of Han Chinese ethnicity, residing in Fujian, China.

According to the “Diagnosis, treatment and prevention of severe acute respiratory syndrome coronavirus 2 infection in children: experts’ consensus statement (Fifth Edition)” published in China, the clinical classification of COVID-19 in children includes five categories: asymptomatic infection, mild, moderate, severe, and critical. Mild cases are mainly characterized by acute upper respiratory tract infection, which may be accompanied by fever, fatigue, headache, and muscle pain. General clinical data were obtained from the medical records system.

MRI acquisition

MRI scans were performed using a 3.0T MR scanner (uMR890, United Imaging, Shanghai, China) equipped with a 64-channel head coil at Xiamen Children’s Hospital. To minimize head motion, foam padding was applied, and noise-canceling headphones were provided to reduce machine noise interference. High-resolution 3D T1-weighted structural MRI data were obtained using a GRE-FSP sequence with the following parameters: repetition time (TR) = 6.7 ms, echo time (TE) = 2.3 ms, flip angle = 8°, acquisition matrix = 320 × 300, field of view (FOV) = 256 × 240 mm, slice thickness = 0.8 mm, and spatial resolution = 0.8 × 0.8 × 0.8 mm. Diffusion-weighted images were acquired using an EPI-DTI sequence across three b-values (500, 1000, and 3000 s/mm²) with 9, 12, and 48 non-collinear directions, respectively. The scan parameters were: TR = 3016 ms, TE = 77.1 ms, flip angle = 90°, acquisition matrix = 140 × 140, FOV = 210 × 210 mm, slice thickness = 1.5 mm, and spatial resolution = 1.5 × 1.5 × 1.5 mm.

Image processing

Structural images were processed using the FreeSurfer package (version 7.4.0; https://surfer.nmr.mgh.harvard.edu/) for cortical reconstruction and volumetric segmentation. The technical details of these procedures have been previously described in prior publications24,25. Cortical metrics, including cortical thickness, surface area, volume, and the Local Gyrification Index (LGI), were measured using surface-based morphometry principles. To enhance the robustness of the analysis, individual cortical metric maps were smoothed with a Gaussian kernel of 5 mm full-width at half maximum. The brain was parcellated into 68 cortical regions using the Desikan-Killiany (DK) atlas, and the mean LGI for each region was calculated26. Given the developmental characteristics of our pediatric sample (ages 1–6), additional steps were taken to ensure the accuracy of segmentation and parcellation results. Beyond the standard FreeSurfer pipeline, we also processed the data using Infant-FreeSurfer for validation purposes. All segmentation results were manually reviewed by two experienced pediatric neuroradiologists, who evaluated the anatomical plausibility of the parcellations. Our findings showed that for the 1-6-year-old cross-age sample, the standard FreeSurfer pipeline produced more stable and accurate results. Consequently, we based subsequent analyses on the segmentation results obtained from the standard FreeSurfer pipeline.

Diffusion images were preprocessed using MRtrix3 (version 3.0.4; https://www.mrtrix.org/) and DSI Studio (version 2022; https://dsi-studio.labsolver.org/) following the procedures outlined below27: (1) quality inspection; (2) conversion of NIFTI to the native MIF format; (3) denoising using random matrix theory28; (4) removal of Gibbs ringing artifacts29; (5) motion and distortion correction30; (6) bias field correction; (7) conversion of MIF back to NIFTI for further processing in DSI Studio; (8) creation of a mask to exclude background regions and enhance reconstruction accuracy; and (9) DTI-based reconstruction to characterize the major diffusion direction of the fibers.

Construction of structural covariance networks

SCNs based on the LGI were constructed using the Graph Analysis Toolbox (GAT)31,32. For each group, a 68 × 68 correlation matrix was established based on the DK atlas, with connection strength defined as Pearson correlation coefficients between LGI histograms of each pair of brain regions, adjusted for age, gender, and estimated Total Intracranial Volume (eTIV) (Fig. 1A, B). The SCNs were subsequently binarized using density thresholds ranging from 0.24 to 0.38 (with an increment of 0.01), ensuring that both HCs and COVID networks contained the same number of nodes and edges at each density (Fig. 1C, D). The minimum density (Dmin) was set at 0.24 to maintain adequate connectivity and avoid network fragmentation in both groups. The density ranged from Dmin to the maximum density at which the networks in both groups exhibited a small-world index (σ) greater than 1.2.

Fig. 1
figure 1

Correlation and binary adjacency matrices for healthy controls (HCs) and COVID-19 patients. (A) Correlation matrix for HCs, (B) correlation matrix for COVID-19 patients, (C) binary adjacency matrix thresholded at Dmin (0.24) for HCs, and (D) binary adjacency matrix thresholded at Dmin (0.24) for COVID-19 patients. The color bar represents the correlation coefficient, indicating the strength of the connections between brain regions.

Full size image

Global and nodal network metrics were calculated in accordance with established definitions from previous studies33,34. The global metrics included the clustering coefficient (Cp), characteristic path length (Lp), normalized clustering coefficient (γ), normalized characteristic path length (λ), σ, global efficiency (Eg), transitivity, and modularity. Nodal metrics included the nodal clustering coefficient, nodal degree, nodal betweenness centrality, and nodal local efficiency. Network hubs were identified as regions where the nodal betweenness centrality was at least 2 standard deviations larger than the average value31. The hubs of both groups were visualized using the BrainNet Viewer35.

Calculation of DTI-ALPS index

Fractional anisotropy color maps were generated using DSI Studio. At the level of the lateral ventricle body, regions of interest (ROIs) with a 3 mm diameter were placed on bilateral projection and association fibers by an experienced neuroradiologist. Fiber orientation and diffusivities along the x-axis (right-left; Dxx), y-axis (antero-posterior; Dyy), and z-axis (infero-superior; Dzz) within the ROIs were obtained. The DTI-ALPS index for the left hemisphere (left ALPS) and right hemisphere (right ALPS) was calculated separately using the formula: DTI-ALPS index = mean (Dxxproj, Dxxassoc)/ mean (Dyyproj, Dzzassoc)36. A higher DTI-ALPS index indicates better glymphatic function.

Estimation of CP volume

Bilateral CP volumes were obtained from the automatic segmentation of 3D T1-weighted images using FreeSurfer. To minimize interindividual variability, CP volume was normalized to the eTIV and multiplied by 1000, as recommended in prior studies37. A higher CP volume is indicative of increased CSF production.

Statistical analyses

Statistical analyses were conducted using IBM SPSS 23.0. Continuous variables are presented as mean ± standard deviation or median (IQR), while categorical variables are expressed as numbers (percentages). Group comparisons of continuous variables were performed using the independent samples t-test for normally distributed data, Welch’s t-test for normally distributed data with unequal variances, or the Mann-Whitney U-test for non-normally distributed data. Categorical variables were compared using Chi-square tests. Statistical significance was determined using a two-tailed P-value of < 0.05.

Surface-based cortical metrics were compared between groups using vertex-wise general linear models in FreeSurfer, adjusted for age, gender, and eTIV. Multiple comparisons were controlled using vertex-wise P < 0.01 and cluster-wise P < 0.05. The area under the curve (AUC) for network parameters across all densities was calculated and compared between groups. A non-parametric permutation test with 1,000 repetitions was applied to assess significant differences in global and nodal network metrics. Comparisons of network metric were performed using the GAT toolbox, with statistical significance defined as P < 0.05, corrected for false discovery rate (FDR).

Results

Demographics

A total of 19 children diagnosed with COVID-19 and 22 HCs were included in the study. No significant differences were observed between the groups in terms of age (Z = -0.894, P = 0.379) or gender (χ2 = 3.107, P = 0.098). Additionally, there were no significant differences in eTIV between the groups (1,243,079.3 ± 166,049.0 mm3 vs. 1,183,346.8 ± 161,198.2 mm3; t = 1.164, P = 0.251). Among the COVID-19 patients, 14 (73.7%) experienced fever, 10 (52.6%) had a cough, 2 (10.5%) reported a headache, 2 (10.5%) had diarrhea, and 1 (5.3%) experienced myalgia.

Between-group comparisons of cortical metrics

Vertex-wise comparisons of cortical metrics are summarized in Table 1. Compared to the control group, the COVID-19 patients showed a significant increase in cortical area, volume, and LGI in the left superior parietal cortex. Additionally, cortical thickness was significantly increased in the left lateral occipital cortex (Monte Carlo correction, vertex-wise cluster threshold = 0.01, cluster-wise P-threshold = 0.05).

Table 1 Regional changes in cortical metrics in the children with COVID-19 compared with HCs.
Full size table

Between-group differences in global and nodal network parameters

To further investigate group differences, we assessed global and nodal network parameters based on cortical LGI. Changes in global network parameters across densities (ranging from 0.24 to 0.38) for the COVID and HCs groups are illustrated in Fig. 2. The COVID group exhibited a significantly lower Cp compared to the HCs, with significant difference observed at several density thresholds (P < 0.05, FDR-corrected; Fig. 2B). No significant differences were observed between the groups for the remaining parameters (Fig. 2). Notably, both groups’ SCNs displayed small-world properties (σ > 1.2; Fig. 2I). The AUC comparisons further confirmed a significantly lower clustering coefficient in the COVID group (P = 0.028, FDR-corrected), while no other measures reached statistically significant between the groups.

Fig. 2
figure 2figure 2

Global network parameters for HCs and children with COVID-19, including between-group differences. (A,B) Clustering coefficient, (C,D) characteristic path length, (E,F) normalized clustering coefficient, (G,H) normalized path length, (I,J) small-world index, (K,L) global efficiency, (M,N) transitivity, and (O,P) modularity of the HCs and COVID-19 networks. Red * markers outside the 95% confidence interval (dashed lines) denote network densities with significant difference (P < 0.05). Positive values indicate higher densities in children with COVID-19 compared to HCs, while negative values indicate the opposite.

Full size image

AUC analyses of four nodal network measures—clustering coefficient, degree, betweenness centrality, and local efficiency—were also performed (Fig. 3). However, no significant differences between the groups were found after FDR correction.

Fig. 3
figure 3

Between-group differences in network parameters across a range of network densities: (A) normalized clustering coefficient, (B) normalized degree, (C) normalized betweenness, and (D) normalized local efficiency. Red * markers outside the 95% confidence interval indicate regions where the difference between the two groups at a given density is statistically significant.

Full size image

Network hubs

Group-specific network hubs for the HCs and COVID groups are shown in Fig. 4. In the HCs, network hubs were identified in the left precuneus, left rostral middle frontal cortex, right pericalcarine cortex, right superior parietal cortex, and right insula. In contrast, the number of hubs in the COVID group was reduced, with hubs identified only in the right postcentral and right precentral regions.

Fig. 4
figure 4

Network hubs identified in the HCs (A) and COVID group (B). Five network hubs are highlighted in red for the HCs, while two network hubs are highlighted in blue for the COVID group.

Full size image

Between-group comparisons of DTI-ALPS index and CP volume

Figure 5; Table 2 illustrate the differences in bilateral DTI-ALPS index and CP volume between the COVID-19 patients and HCs. Compared to HCs, the COVID group showed significantly higher CP volumes in the left hemisphere (0.36 ± 0.13 vs. 0.28 ± 0.07, P = 0.041) and in total CP volume (0.40 ± 0.13 vs. 0.33 ± 0.07, P = 0.043). However, no significant differences were observed between the groups in the left ALPS index (1.39 ± 0.17 vs. 1.31 ± 0.12, P = 0.081) or in right CP volume (0.37 ± 0.09 vs. 0.44 ± 0.13, P = 0.070).

Fig. 5
figure 5

Differences in DTI-ALPS index and choroid plexus (CP) volume between healthy controls (HCs) and children with COVID-19.

Full size image
Table 2 ALPS index and choroid plexus (CP) volume.
Full size table

Discussion

To the best of our knowledge, this study is the first to investigate the alterations in cortical morphometric parameters, LGI-based geometric networks, and glymphatic function in children with COVID-19. Our main findings include: (1) Children with COVID-19 exhibited significantly larger cortical area, volume, and LGI in the left superior parietal cortex, along with increased cortical thickness in the left lateral occipital cortex; (2) A reduction in the Cp of global network properties and alterations in network hubs were observed in young children with COVID-19, despite the absence of significant differences in regional network parameters; and (3) Both total and left CP volumes were higher in children with COVID-19 compared to HCs at least two weeks after recovery. These results suggest that abnormal cortical metrics, altered geometric networks, and potential glymphatic dysfunction may be consequences of COVID-19 in the developing brain.

The precise mechanisms through which COVID-19 induces neurological damage remain unclear. Current theories propose that SARS-CoV-2 may cause brain injury through direct invasion of the nervous system, disruption of the blood-brain barrier, vascular endothelial damage, and autoimmune responses38. The resulting dysregulated immune activation, neuroinflammation, and cytokine storms are believed to contribute to neurological dysfunctions in pediatric patients, with the severity of these effects correlating with the intensity of the disease39,40,41. These insights are consistent with our findings and provide a possible explanation for the observed alterations in cortical morphometric parameters and glymphatic function in children with COVID-19.

The observed increases in cortical area, volume, and LGI in the left superior parietal cortex may indicate heightened neuroplastic responses to COVID-19 in young children, potentially serving as compensatory mechanisms during critical developmental stages when this region plays a pivotal role in spatial cognition and sensorimotor integration42. Similarly, cortical thickening in the left lateral occipital cortex could reflect an adaptive neuroprotective response to preserve visual processing capabilities during the early stages of recovery43. These neurostructural patterns starkly contrast with established findings in adults, which demonstrate COVID-19-associated cortical thinning and gray matter reduction across multiple cohorts7,10,11,13,44. Systematic comparison with adult literature reveals partial overlap in the brain regions involved, despite the opposing morphological changes. Notably, studies have documented COVID-19-related cortical thinning in the left lateral occipital cortex among young adults with mild infections, while gray matter volume reductions in the left superior parietal cortex are correlated with cognitive deficits in adult patients11,45. However, no significant alterations in the cortical area of the left superior parietal cortex have been observed in adult patients11,46. Two key factors may explain these discrepancies. First, our study focused on a pediatric population aged 1–6 years, during a developmental stage when the brain may be more responsive or capable of compensating for the effects of COVID-19. Second, the cortical changes observed in our study represent an early or transient phase of the brain’s response to infection, differing from the long-term effects reported in the adult literature. Longitudinal studies are necessary to determine whether these cortical changes represent protective adaptations, pathological markers, or precursors to later neurodevelopmental consequences. Such studies would provide a clearer understanding of the long-term impact of COVID-19 on pediatric brain development.

Notably, the abnormal brain regions identified in our study are predominantly localized to the left hemisphere. While the neurobiological mechanisms underlying this lateralization pattern require further investigation, our findings align with recent neuroimaging evidence on COVID-19. Specifically, Douaud et al. reported left-lateralized alterations in cortical thickness in adults with post-COVID conditions7. Similarly, Arrigoni et al. observed leftward-asymmetric microstructural abnormalities in structural connectivity among COVID-19 patients with persistent olfactory or cognitive deficits47. Furthermore, this lateralization pattern may not be limited to adults, as demonstrated by Invernizzi et al., who identified left-lateralized functional and structural brain alterations in adolescents following SARS-CoV-2 infection13. The consistency of these findings across diverse patient populations and imaging modalities suggests that SARS-CoV-2 infection may disproportionately affect left-hemispheric circuits involved in multisensory integration and higher-order cognitive processes.

In recent years, brain connectomics has emerged as a crucial tool for studying brain pathology22,48. To explore the neuroanatomical changes associated with mild COVID-19, we analyzed SCNs based on cortical LGI, revealing potential alterations in cortical connectivity linked to the disease. Our findings indicate that the COVID group exhibited a significantly lower Cp in global network properties compared to HCs, suggesting disruptions in local efficiency and possible segregation of cortical regions following COVID-19 (Fig. 2A, B). This may be attributed to neuroinflammation or other post-infectious processes49. Similar patterns have been observed in SCN studies of patients with ischemic moyamoya disease and individuals with autism spectrum disorder, highlighting the relevance of these findings50,51. Despite these disruptions, small-world properties were preserved in both groups, consistent with the brain’s inherent resilience. Previous research on neurological conditions has noted that global efficiency often remains intact despite localized disruptions, suggesting the brain’s capacity to maintain overall connectivity52. However, we did not observe significant differences in nodal network properties (Fig. 3), indicating that the effects of COVID-19 on brain connectivity may be subtle or dispersed across multiple nodes. This could contribute to the observed global network alterations without significantly affecting individual node characteristics. It is also possible that these nodal changes were too small to detect given the current sample size or statistical power. Additionally, the study identified network hubs in HCs distributed across the left precuneus, left rostral middle frontal cortex, right pericalcarine cortex, right superior parietal cortex, and right insula regions—areas typically involved in sensory integration and complex cognitive functions53,54,55,56. In contrast, the reduced hubs in the COVID group, concentrated in primary motor and sensory areas, suggest that COVID-19 may impair higher-order cognitive processing and multimodal integration57,58.

The glymphatic system, which plays a crucial role in clearing brain waste through CSF circulation, is essential for maintaining the health of the nervous system59,60. A lower CP volume and a higher DTI-ALPS index have been associated with better glymphatic function61,62,63. In our study, we found significantly higher CP volumes in the left and total brain regions of COVID-19 patients compared to HCs, along with a trend towards a reduced left DTI-ALPS index in the patient group. The increased CP volume in the COVID group may indicate an underlying inflammatory response or altered CSF dynamics, potentially acting as a compensatory mechanism or as a consequence of viral invasion and immune activation64,65. Although the reduction in the left DTI-ALPS index was not statistically significant, it may suggest early-stage glymphatic dysfunction. Previous studies have reported significantly lower bilateral DTI-ALPS indices in COVID-19 patients four months post-infection compared to HCs12. The discrepancies in our findings could be attributed to several factors: first, our study may have captured an earlier stage of recovery; second, our cohort had fewer neurological symptoms and milder disease severity, which could influence glymphatic function; and third, limitations in sample size and imaging techniques may have affected the detection of statistically significant differences.

This study has several limitations. First, the small sample size and cross-sectional design may impact the robustness of our results and limit our ability to explore the long-term effects of COVID-19 on the developing brain. Future studies with larger sample sizes and longitudinal follow-up are needed to provide deeper insights. Second, the SCNs were constructed as group-level networks rather than individual networks for each participant. In future work, we plan to create individual-level networks, such as morphological similarity networks and morphometric inverse divergence networks, to provide a more detailed and personalized representation of network characteristics66,67. Third, using a more fine-grained atlas specifically designed for infants may offer more detailed information about the SCNs in children with COVID-19. Fourth, we did not collect clinical data from patients, which prevented us from analyzing correlations between clinical variables and cortical metrics, topological properties, and glymphatic system-related metrics. Lastly, as this is an exploratory study, we did not adjust for covariates such as age and gender in the group comparisons of glymphatic system-related metrics, which is another potential limitation.

In conclusion, our study provides the first report of changes in cortical metrics in young children with COVID-19, despite no observed differences in regional SCN parameters. We also employed a non-invasive multiparametric approach to assess the glymphatic function in COVID-19 patients. We speculate that viral-triggered neuroinflammation and immune response contribute to neurotoxic damage in brain regions involved in cognitive processing, ultimately resulting in alterations in cortical morphology and glymphatic function.