Diffusion Kurtosis Imaging of neonatal Spinal Cord in clinical routine

2.1 Subjects

Infants whose data have been used to disclose each step of the pipeline have been enrolled since August 2019 and scanned with 3.0 T MR scanner using a 32-channel head array coil (Ingenia Cx, Philips, Best, the Netherlands) at the Neuroradiology Unit of Giannina Gaslini Children’s Hospital of Genova. Conventional MRI and DKI were performed in 17 pre-term infants (28.1 to 36.7 weeks Gestational Age (GA); scanned at Term-Equivalent Age (TEA)). Details about subjects’ demographics are reported in Table 1.

This single-center study was carried out in accordance with the recommendations of “Comitato Etico Regione Liguria, Genova, Italy” with written informed parental consent obtained for each infant prior to examination in accordance with the Declaration of Helsinki. Subjects were spontaneously breathing during examination; free-flowing oxygen was administered for all the duration of MRI session if necessary. Throughout the course of the examination, newborns were subjected to constant monitoring of the oxygen saturation and heart rate, by pulse oximeter and three-electrode electrocardiographic monitoring, respectively.

Exclusion criteria included obvious motion artifacts, oblique positioning, an incomplete imaging process or a low SNR. We opted for focusing exclusively on upper SC portion, by only acquiring the cervical district (C1-C7), affected by a significantly lower degree of physiologic motion than the dorsal and lumbar portions (38).

In consensus with a board-certified pediatric neuroradiologist, we performed QC for each of the pipeline’s steps.

2.2. Full pipeline description

Our pipeline integrates MRtrix3 (v.3.0.1) (39) for setting of dMRI acquisition sequence, Spinal Cord Toolbox (SCT, v. 5.3.0, https://github.com/neuropoly/spinalcordtoolbox) (40) for all processing steps specific to the SC, and Diffusion Imaging in Python (Dipy, v.1.4.0) (41) for denoising as well as computation of diffusion metrics.

Output of key processes, such as motion correction, segmentation and registration with atlas, can be checked through an interactive SCT Quality Control (QC) module, which automatically generates reports consisting in HTML files, containing a table of entries and allowing to show, for each entry, animated images (background with overlay on and off) for data quality validation.

In our methodological pipeline we have opted for mainly relying on SCT, being currently the only existing fully-comprehensive, free and open-source software dedicated to the processing and analysis of multi-parametric MRI of the spinal cord (42–53) successfully employed in a plethora of clinical applications concerning adult SC (54–63).

An overview of our image processing pipeline highlighting key features is shown in Figure 1. Since SCT algorithms are validated in adult imaging, we specifically customized each processing step to our neonatal scans. Our pipeline thus represents, to the best of our knowledge, the first semi-automated ad-hoc procedure for imaging of neonatal spine. A fully automatic workflow is not feasible here: acquisition time constraints, available scanner features, and subsequent image quality require inevitable although minimal and highly reproducible manual interventions.

• Download figure
• Open in new tab

Figure 1. Overall Processing Pipeline:

designed pipeline allows complete handling of DKI scan of neonatal Spinal Cord from acquisition setup to preprocessing, processing and postprocessing

2.3. Customized acquisition setting

In order to minimize macroscopic movement artifacts, all recommended guidelines for pediatric imaging have been adopted. So as to protect infants from acoustic disturbances caused by MR sequences, we resorted to baby earmuffs and silicone paste for hearing aids. Furthermore, we avoided most of motion by swaddling infants and by placing airbags around the baby’s head. In addition, protective pads have been placed between the magnet and the patient. All these contribute to create a comfortable and warm rest environment, minimizing the chance of free movements.

MRI was performed when possible during spontaneous sleep by exploiting the administration of breast milk or formula about 30 minutes before the start of the exam. In case of spontaneous sleep failure, in order to minimize macroscopic movement artifacts, the instrumental examination was performed under mild sedation by orally administering Midazolam at 0.1-0.2 mg/kg diluted in glucose solution 33%, subject to signature of informed consent from parents and applied by expert trained nurses.

Given the lack of a specific acquisition protocol for DKI of neonatal SC, we designed the diffusion-weighting scheme in collaboration with the neuroradiologists at Giannina Gaslini Hospital. One constraint we had to deal with was the impossibility to perform optimized variants of Spin-Echo Echo Planar Imaging (SE-EPI) sequence (i.e. reduced FOV or spatially selective techniques) (20) on Philips Ingenia scanner. Therefore, minimization of scan duration was our main focus in order to suppress motion and fast CSF pulsation artifacts typical of newborns.

We thus tested different versions of diffusion-weighted gradient scheme, adopting optimal trade-off between Fiber Orientation Distributions profile (FODs, estimated with Mrtrix3 using Multi-Shell Multi-Tissue Constrained Spherical Deconvolution (MS-MT CSD)), image quality and scan time.

We generated each multi-shell diffusion gradient table through Mrtrix3 script gen_scheme, taking as inputs the number of phase-encoding directions to be included in the scheme (for most scanners, including ours, typically 1), the b-value of the shell, and the number of directions to include in the shell. This procedure ensures uniform sampling by maximizing uniformity within shells using a bipolar electrostatic repulsion model for optimal angular coverage.

For further reducing acquisition time without significantly affecting image quality, we applied the MultiBand slice acceleration technique (64) (https://www.usa.philips.com/healthcare/resources/landing/compressed-sensecombined).

Moreover, spatial saturation bands have been used to suppress strong signals originating from the subcutaneous lipids to enable volumetric region coverage.

The final version of diffusion acquisition scheme is displayed in Figure 2 as well as reported in Table 2, and includes 6 b=0; 13 b=700 and 13 b=2100 s/mm² for a duration of 4 minutes 30 seconds. This allowed acquisition of high in-plane resolution axial diffusion weighted images, where b=0 scans could be well discriminated from non b=0 volumes and anatomical SC features are sharp.

• Download figure
• Open in new tab

Figure 2. 3D view of final diffusion acquisition scheme:

directions of diffusion-sensitizing gradients relative to each b-value are displayed in three different colors as reported in the legend. Units are in s/mm². Markers indicates polarity: dots are the polarities in the set, asterisks their opposite Along with dMRI, we also acquired a high-resolution structural image as anatomical reference. The definitive MRI protocol thus consisted in a Turbo Spin Echo (TSE) 3D T1-weighted image followed by a DKI series whose details are listed in Table 1.

This is the best possible solution that combines conflicting requirements of complex multi-shell acquisition and minimized scan time.

2.4. Preprocessing

2.4.1. Denoising

SC imaging is characterized by low SNR, which can hamper accurate, repeatable, quantitative measurements. Moreover, models such as DKI are susceptible to noise and signal fluctuations often leading to degeneracies in estimation of derived parameters. SNR further lowers in case of neonates, due to relatively high overall free water content, and denoising approaches based on Principal Component Analysis (PCA) are inapplicable due to a reduced number of diffusion gradient directions.

Therefore, we adopted Patch2Self, a recently proposed self-supervised learning denoising method that outperforms existing non supervised methods. (65).

A unique advantage of Patch2Self is the lack of requirement for selecting or calibrating an explicit model either for noise or diffusion signal, so that it can be applied at any step in the pre-processing pipeline. The only assumption it relies on is randomness and uncorrelation of noise across different gradient directions. Its framework consists in holding out one volume and using patches from all other volumes to predict the center of the patches of the held-out volume using a regressor. This denoiser has already showed a significant improvement in repeatability and conspicuity of pathology in diffusion volumes and quantitative DTI metrics for adult SC (66).

Here, we chose to apply Patch2Self as first preprocessing step on raw data since it showed to offer highest SNR. The method is implemented in Dipy v.1.4.0 and applied with Ordinary Least Squares (OLS) regressor, since recommended for SC imaging (Figure 3).

• Download figure
• Open in new tab

Figure 3. Visual inspection of denoising:

The denoising of Patch2Self is compared against the original noisy image along with their corresponding residuals. Notice that Patch2Self does not show any anatomical structure in the corresponding residual plots likely neither introducing structural artifacts

2.4.2. Cropping

SC scans usually include also cerebral areas such as medulla and cerebellum due to their proximity with cervical SC. In order to exclusively focus on the area of interest excluding undesired voxels, as a first preprocessing step, we thus recommend to apply to DKI images SCT function sct_crop_image allowing also to fasten subsequent processing. Lower and higher bounds for cropping along the three spatial coordinates can be specified via command line in order to select the same area of interest (i.e. cSC) for all the cohort, considering FOV positioning is consistent across subjects.

Specifically, in the case of our scans, FOV reduction allowed to exclude upper non-spinal areas (i.e. cerebellum) as well as lower spinal levels whose corresponding slices are not usable due to poor image quality (Figure 4a).

• Download figure
• Open in new tab

Figure 4. Preprocessing:

DKI scan through preprocessing steps for one example subject: (a) FOV reduction; (b) Motion Correction; (c) Segmentation: Deep Learning Segmentation algorithm generally achieves satisfactory results in SC detection; (d) Example of artefactual slice due to a poor fat saturation, causing the fat to alias on the spinal cord area and (e) requiring manual correction of segmentation

2.5. Processing

2.5.1. Vertebral Labeling

After segmentation, labeling of vertebral levels or discs is the second mandatory step in order to match the template to the subject’s MRI (template registration).

Two vertebral levels are necessary for registering data to the template. Each of these two landmarks consists of a voxel placed in the middle of the SC, at the level of the corresponding mid-vertebral body, and assigned a relative number starting from 1 for C1 vertebra. However, SCT recently introduced the possibility to alternatively use inter-vertebral disc labels with the analogous procedure of reference numbered voxels.

We perform this step on 3D T1w images in order to achieve better accuracy given their higher overall quality and contrast compared to DKI ones, where vertebral discs are not clearly identifiable.

Labeling from 3D T1w anatomical image is possible as it turned out to match relatively well along the superior-inferior (z) axis, the target direction of disc labeling, with the DKI scan (not along the Anterior-Posterior (AP) or Right-Left (RL) direction, Figure 5).

• Download figure
• Open in new tab

Figure 5. DKI scan overlaid on structural 3dT1w image:

while both images are clearly not registered along the antero-posterior direction due to the very strong susceptibility artefact, the z-location is similar: see how the bottom tip of the cerebellum is consistent for the two scans

Vertebral labeling is typically done using an automatic method sct_label_vertebrae, that finds C2-C3 disc, and then locates neighbouring discs using similarity measure with the PAM50 template at each specific level (70).

Default SCT procedure sct_label_vertebrae fails in automatically detecting C2-C3 vertebral disc once again because of the small size of spines at issue and low image contrast compared to adults.

Therefore, we manually create labels with the command sct_label_utils through interactive viewer option provided by SCT (flag -create-viewer) with little to no waste of time.

Specifically, vertebral labeling was created at the posterior tip of the top of C1 vertebra and at C3-C4 disc, centered in the cord. Manual intervention only took a few seconds per subject (Figure 3S).

2.5.2. Registration to PAM50 atlas

Registration between subject’s diffusion and atlas space is a very demanding task in case of neonatal imaging given the lack of a specific pediatric atlas compatible with SCT (one is currently under creation, https://github.com/neuropoly/spinalcordtoolbox/issues/2530). We thus use PAM50 atlas (71), an adult template for MRI of the full SC and brainstem in the same coordinate system as the ICBM152 (MNI) brain template, allowing to conduct simultaneous brain/spine studies. It consists of a T1w, T2w, T2*w, white and gray matter probabilistic atlas and white matter atlas of tracts as well as probabilistic labeling of spinal levels. The template has been constructed from straightened SC for facilitating registration and visualization of results.

sct_register_to_template is the main command for registering one subject to the template and vice versa, since it outputs the forward and backward warping fields. We choose subject’s native diffusion space as target of registration transforms since the straightening required by opposite strategy would cause through-plane interpolation errors which would bias following extraction of diffusion measures (72).

Moreover, we suggest employing T1w atlas image for its better contrast similarity with DKI scan compared to T2w.

Application of default command does not produce satisfactory results, stressing the need to tweak all input parameters to deal with our particular contrast and resolution. Given the presence of artifacts and some inherent features (e.g., low CSF/cord contrast) that could compromise the registration, we use SC segmentation as input for the algorithm in order to ensure maximum robustness.

Registration is then built through multiple steps by increasing the complexity of the transformation performed in each step (starting with large deformation with low degree of freedom and finishing with local adjustment). Specifically, the first step consists in vertebral alignment, that is vertebral level matching between the subject and the template on the basis of posterior edge of the intervertebral discs provided by previous manual vertebral labeling. Second step is slice-wise center of mass alignment between the two images, using centermass algorithm instead of default centermassrot (which also includes rotation alignment) because the cord is quasi-circular and cord angle estimation is not reliable here. The third step is R-L scaling along x axis followed by A-P alignment to match segmentation borders along y axis, with the ultimate aim of accommodating the very small SC size. Finally, iterative slice-wise non-linear registration is performed through non-linear symmetric normalization regularized with b-splines (73) using information from comparison of Cross Correlation metric (CC) between the two images, which allows SC shape refinement. Once the algorithm completed, one can assess the quality of registration through visual evaluation and inspection of QC module, and thus warp the template and all its objects to each subject’s DKI image (Figure 6).

• Download figure
• Open in new tab

Figure 6. Registration with PAM50 atlas and ROI detection through atlas-based approach:

(a) PAM50 atlas’ cord segmentation binary mask, (b) WM and (c) GM probabilistic masks warped to subject’s DKI motion-corrected mean image

Current selection of parameters and steps successfully worked for our scans since atlas registration algorithm robustly achieved convergence as verified through inspection of QC feature.

2.6 Postprocessing

Thanks to this atlas-based analysis approach, it is possible to perform cord-specific quantification of diffusion metrics through sct_extract_metric command, also restricted to specific Regions of Interest (ROIs; labels used by default are taken from the PAM50 template, e.g. WM tracts, flag -l), vertebral levels (flag -vert) or slice (flag -z), according to the specific clinical needs concerned.

Along with WM and GM probabilistic masks as a whole (Figure 6b,c), normally investigated in medical practice, one can carry out ROI detection also in specific tracts according to the clinical question (fifteen WM tracts and three GM regions available in total for each side).

In our example, neither DKI nor structural images ensured sufficient WM-GM-CSF contrast to perform any manual detection of ROIs in contrast to high-contrast PSIR image of (82), whose acquisition time would be too long for neonates. Therefore, we exploited good registration outcome for automatic delineation of ROIs through atlas-based approach.

We opted for using lateral CSTs as ROIs for consistency with (82) - though grouping together left and right sides in order to gain robustness by increasing volume fraction as suggested in (40) - as well as WM and GM.

We then computed averages of diffusion measures across C1-C4 vertebral levels in order to avoid a subject-based bias: outside of these levels the registration is inaccurate and/or MRI signal may be corrupted. We thus checked through QC module the correctly segmented slices corresponded to the same vertebral levels across subjects, starting from the first slice containing only SC (excluding cerebellum, Figure 2Sc).

Moreover, estimation of DTI and MSDKI weighted average metrics was limited to those slices where spinal cord segmentation is accurate: outside segmentation mask, metrics would indeed be irrelevant. This was obtained by multiplying segmentation mask by specific WM, GM and CSTs atlas labels. We quantified diffusion metrics using Weighted Average (WA) estimation to minimize PVE avoiding bias into resulting metrics by the surrounding tissues (e.g. CSF). This is one of the recommended methods especially in case of noisy images and small tracts as in our case. We assessed associated voxel fraction to quantify the reliability of our diffusion measures: as demonstrated in (40), having at least 240 voxels results in an error smaller than 1%, while having 30 voxels results in an error inferior than 2%. In this example, the metrics were computed based on average 178.3, 50.5 and 31.5 voxels in WM, GM and CSTs respectively, thus assuring sufficient accuracy of estimates.

2.7 Case Study

Periventricular WM Injury (PMWI) is the most frequent type of brain lesion in preterm infants, and the spatial extent and location of WM injury correlate with distinct clinical outcomes, including cerebral palsy and motor impairment (83).

Given the strong association of WM injury with the motor function development of preterm neonates, we assessed the hypothesis that the presence of periventricular punctate WM lesions identified on MRI at TEA could be associated with regionally specific alterations in cSC microstructure.

A similar approach was already used by (82) to characterize cSC microstructural abnormalities in a cohort of adult patients with previous unilateral ischemic stroke in the vascular territory of the middle cerebral artery (MCA). DTI and DKI diffusion measures in cSC resulted to be a valuable imaging marker for predicting clinical outcome. In particular, significant reduction of FA and MK was observed in the affected lateral WM bundle of the cSC, correlating with the severity of motor dysfunction.

Accordingly, the ultimate goal of our study was to verify whether the presence of periventricular WM lesions affects the cSC tracts development. Specifically, we aimed to compare DTI and MSDKI measures of cSC in two groups of preterm neonates: (i) with punctate Periventricular White Matter lesions (PWMI), and (ii) with normal brain MRI (controls).

3.1. Population size and classification

In order to investigate clinical differences among acquired subjects, we grouped infants as follows: i) 9 subjects with punctate PWMI and ii) 9 subjects with normal brain MRI, used as control group. At QC phase, in accordance with the expert neuroradiologist, we opted for excluding one control subject due to excessively poor image quality (i.e. signal leakage at C1-C3 level, Figure 2Sa,b). Therefore, the final number of subjects under analysis amounted to 9 and 8 infants for the two subsets, respectively.

3.2. The role of denoising

As mentioned above, neonatal imaging is inherently affected by low SNR and sensitive to imaging artifacts. Proper denoising of scans is therefore a crucial step in the processing pipeline. Above all, we thus focused on quantitatively assessing the contribution of Patch2Self denoiser on subsequent analysis.

Firstly, we computed average SNR on b=0 images, namely thermal noise due to the diffusion-related signal attenuation, across all slices belonging to C1-C4 district of our interest.

As expected, we found a remarkable increase in mean SNR after applying Patch2Self (5.88 ± 1.41 vs 14.64 ± 4.53) as displayed in Figure 7.

• Download figure
• Open in new tab

Figure 7. Effects of Patch2Self denoising on thermal noise.

Average Signal-to-Noise Ratio (SNR) computed on b=0 images increases, across C1-C4 vertebral levels under analysis, when including denoising with Patch2Self algorithm in the processing pipeline

We then inspected impact of denoising on microstructure model fitting, a critical step often leading to degenerate parameter estimates due to the low SNR of dMRI acquisitions.

Specifically, we applied DTI and MSDKI model on raw and denoised data, resorting both to traditional Marchenko-Pastur PCA (MP-PCA) and to Patch2Self method. We opted for comparing our denoising procedure with MP-PCA since it represents the current state-of-the-art unsupervised method for denoising DWI (84). In order to compare the goodness of each fit, we performed a k-fold cross-validation (85) across an axial slice (1312 voxels) of masked data for an example subject pulled out from our dataset. In Figure 8, we depict the improvement of the metric by simply subtracting the goodness-of-fit scores of fitting noisy data from Marchenko-Pastur and Patch2Self denoised data for both DTI and MSDKI models. Patch2Self performance is not significantly better than MP-PCA in the goodness-of-fit for DTI model (two-sided t-test with Bonferroni correction, p=1). On the contrary, this state-of-the-art denoiser shows a significant improvement for MSDKI (p=1.825e-04), proving once more to be particularly suitable for DKI model.

• Download figure
• Open in new tab

Figure 8. Violin plots quantifying the increase in metric after fitting downstream DTI and MSDKI models for one axial slice in one example subject.

The improvements in each case are plotted by subtracting the scores of model fitting on noisy data (Raw) from of fitting each denoised output. Note that the consistency of microstructure model fitting on Patch2Self (P2S) denoised data is higher than that obtained from Marchenko-Pastur (MPPCA), especially as regards MSDKI model. ns stands for non significant, *** stands for 1.00e-04 < p <= 1.00e-03 in two-sided t-test with Bonferroni correction

3.3. Trend of diffusion measures

For consistency with (82), we considered above all trend of cSC FA and MSK measures in the two groups of neonates as displayed in Figure 9. For the sake of completeness, in Figure 4S we also include remaining DTI metrics.

• Download figure
• Open in new tab

Figure 9. Extraction of diffusion measures within specific ROIs:

(a) White Matter (WM), (b) Gray Matter (GM) and (c) Cortico-Spinal Tracts (CSTs) ROIs overlaid on DKI motion corrected image; (d) Scatter plots of FA and MSK in group subjects across aforementioned ROIs: colored spots indicate single subject’s value for each metric; as reported in the legend, controls’ measures are in blue, whereas Periventricular White Matter Injury (PWMI) group’s in red. Units for MSK are in mm²/s, while FA is dimensionless. Error bars displaying mean (diamond) and standard deviation (bars) are overlaid on scatter plots.

To an initial evaluation based on the limited sample size available, considering the variation of average metrics between controls and patients in WM, we can notice an increase in MD, AD, RD, parallel to an overall decrease in FA and MSK in preterm neonates with PWMI.

Furthermore, range of values of DTI measures are consistent with normative values on healthy pediatric SC (30). MD, AD and RD values are higher while FA values are lower compared to equivalent measures on older cohorts of patients (i.e. children/ adolescents, (22,23,25,26)). This may be partially due to sensitivity of FA to denoising, which can imply a reduction in this metric. Moreover, this trend is in line with the simultaneous age-related decrease in MD, AD and RD and increase in FA metrics reflecting progressive maturation, myelination and fiber packing and thickening within the SC, similar to that observed in the brain (86,87).

Conversely, definition of a normative variation of DKI measures across age from newborns to adults will be feasible after further investigations from early stages of development.

4.1. Research question

In the present study we tested the pipeline of SC DKI analysis in a group of neonates with PWMI, a form of mild WM injury frequently diagnosed in preterm infants. PWMI are seen at brain MRI as small, focal, multiple alterations of signal intensity (high on T1 and/or low on T2) in periventricular WM. Long-term outcome of neonates with PWMI seems to be related to the number of lesions, their pattern, and their localization. Of note, several studies have showed that greater lesion load of PWM and the involvement of frontal WM are associated with higher risk of adverse neurodevelopmental outcome, affecting both motor and cognitive functions (88). Moreover, periventricular WM lesions in preterm neonates are associated with region-specific changes in MD, FA, RD, and AD in several cerebral WM tracts that might explain the abnormal development of long-term neurological functions (89). Specifically, the involvement of pyramidal tract fibers in the periventricular WM has been demonstrated to be a relevant factor for motor dysfunction in children with PWMI (90). In our study, we found that microstructural changes can be detected by using an advanced DKI analysis also in the GM and WM of cSC of preterm neonates with PWMI studied at term-equivalent age, thus demonstrating that DKI parameters could be used as markers to unravel underlying microstructural lesions not visible on conventional structural MRI. Moreover, our preliminary findings confirm the hypothesis that in preterm neonates with PWMI WM microstructure alterations extend beyond the immediate area of periventricular injury, widening distally also in the cSC (91,92). Further analyses on a wider cohort of neonates are necessary to confirm these preliminary results, and specifically to prove if microstructural changes in cSC in preterm neonates with PWMI correlate with long-term neurological outcome.

4.2. Study significance

Here, we present the first application of DKI to neonatal SC through a pipeline able to perform complete processing on a subject within a clinically acceptable time (10 minutes average with the current setup). As regards acquisition setting, we were able to perform a time-consuming technique like DKI using a short diffusion sequence which minimizes patient’s physiological motion and which likely reflects a standard clinical scenario devoid of latest technologies in terms of acquisition sequence optimization. With regard to image processing, we opted for creating this pipeline using SCT since it represents the only existing comprehensive, free and open-source software dedicated to the processing and analysis of SC multi-parametric MRI data. Adaptation of each image processing tools already in use for adult subjects through appropriate tuning of parameters turned out to be feasible. Indeed, it allowed to successfully overcome all the issues mentioned in the ‘Introduction’ section inherent to imaging of SC and exacerbated in case of neonatal setting, even for the most challenging steps like segmentation or registration to atlas. We were thus able to quantify diffusion measures within specific ROIs using an atlas-based approach which presents undisputed advantages compared to usual manual drawing of ROIs. It is automatic and thus highly reproducible, it is not biased by the user experience and knowledge of the anatomy, it is much faster than long and tedious manual delineation of ROIs and it allows to account for PVE. Preliminary application of our method to a limited number of subjects has yielded physiologically plausible findings despite suboptimal acquisition setting. In particular, differences in DTI and MSDKI measures observed between PWMI neonates and controls confirm the hypothesis that PWM injuries in the premature brain can be associated with microstructural alterations of both GM and WM of neonatal cSC at term-equivalent age. These findings also highlight the importance of complementary analysis of DTI and DKI metrics for a more accurate characterization of biological tissues.

4.3 Added value of DKI

Results about feasibility of DTI and MSDKI analysis in neonatal SC with subjects collected so far are preliminary but promising and demonstrate the clinical utility of combining DTI and DKI in the characterization of spinal cord pathologies.

FA reduction parallel to MD increase in patients is an expected finding consistent with existing literature and attributable to degeneration of diffusion barrier and loss of diffusion directionality.

Regarding DKI, it is interesting to note how MSK metric, although yet underused in clinical studies, proves to be more sensitive than FA to microstructural changes in PWMI subjects by further reducing its value in all considered ROIs.

This DKI estimator thus confirms DKI increased sensitivity in capturing alterations related to pathology, also far from the lesion site.

Another crucial food for thought concerns observing how the presence of a WM lesion in the brain causes subsequent alterations not only in cSC WM but also in GM, as evidence of the strong association between brain and spine. In this respect, resorting to DKI measures becomes of utmost importance given kurtosis sensitivity to structural changes in isotropic tissues such as GM.

Indeed, range of variation of MSDKI metric from controls to PWMI is overall higher in GM than that of corresponding DTI measures.

This is an expected finding: DTI is extremely suitable for capturing highly directional diffusion typical of WM bundles but shows low sensitivity to isotropic diffusion, whereas DKI allows specific assessment of non-gaussian diffusion typical of GM areas (93).

Such findings once again stress the importance of combining DTI and DKI metrics as complementary sensitive biomarkers in order to fully exploit the potential of dMRI compared to conventional MRI.

However, a more comprehensive corroboration and explanation of our results is expected after collecting an adequate number of subjects to carry out a robust statistical survey. An in-depth interpretation of the single metrics is out of the scope of this paper. Here, we just dwell on exploring the comparison with the work on adults which served as a starting point.

4.4 Comparison with adult study

Since MSK has proved to be an exact approximation of MK, here we assimilate it with standard DKI metric used in (82).

Both FA and MSK trend agree with just mentioned study. Specifically, they both exhibit a reduction in WM, GM and CSTs in case of pathology (Figure 9).

Since FA is known to be an index of structural integrity (94) and MSK a marker of tissue microstructure complexity (2), obtained findings suggest that, in case of overlying WM brain lesion, a loss of integrity and complexity is registered also in SC WM tracts below with a more isotropic diffusion pattern due to disruption of WM tracts. Also MD, AD and RD follow the same trend, with an increase in case of lesioned subjects, as in (82). This biologically plausible hypothesis already verified for adults would thus also subsist in infants.

4.5. Study limitations

Present pipeline has been designed with current instruments at disposal and with given data whose resolution is high in axial plane but whose contrast is very low due to short acquisition time. As a result, any improvement in acquisition setup our pipeline will bring to even stronger and more comprehensive results.

Major flaws of this procedure consist in basing on an adult atlas, where the exact location of tracts may not perfectly correspond to neonatal images despite the good adjustment of registration parameters. Current pipeline will definitely benefit from introduction of a pediatric atlas into SCT. Moreover, image quality could be further improved: scans we acquired are extracted from routine clinical protocol and consequently prone to noise and artifacts due to short acquisition time dictated by clinical needs and to the lack of specific spatially-selective MR sequences.

Starting from more advanced hardware tools may significantly increase image quality and thus accuracy of estimated metrics. Resorting to optimized acquisition sequences would also allow to increase resolution of HARDI acquisition scheme and thus to exploit all standard DKI measures, which can in turn increase the amount of diagnostic information.

We acknowledge the protocol in use to be on the edge for HARDI schemes required by DKI. However, this represents an attempt to customize advanced dMRI acquisition setting within a clinical routine protocol, already long in itself since made up of multiple MRI sequences in order to increase diagnostic possibilities. Nevertheless, we appropriately addressed this issue at DKI tensor and measures computation phase to ensure reliability and accuracy in their estimates.

4.6. Future developments

Validation of current pipeline can be made by testing it to a larger cohort of subjects, possibly investigating also lower SC tracts, including thoracic and lumbar districts, and extending studies to different clinical cases, preferably focusing on a determined pathology. For example, it would be interesting to explore long-term correlations between DKI measures and specific clinical scores as done in (82), where diffusion measures have been related to motor performance indexes.

A further step may be adapting this analysis pipeline to other promising higher-order diffusion models requiring multi-shell acquisition such as NODDI (95).