Traditionally, research unraveling seasonal neuroplasticity in songbirds has focused on the male song control system and testosterone. We longitudinally monitored the song and neuroplasticity in male and female starlings during multiple photoperiods using Diffusion Tensor and Fixel-Based techniques. These exploratory data-driven whole-brain methods resulted in a population-based tractogram uncovering microstructural sexual dimorphisms in the song control system and beyond. Male brains showed microstructural hemispheric asymmetries, whereas females had higher interhemispheric connectivity, which could not be attributed to brain size differences. Only females with large brains sing but differ from males in their song behavior by showing involvement of the hippocampus. Both sexes experienced multisensory neuroplasticity in the song control, auditory and visual system, and the cerebellum, mainly during the photosensitive period. This period with low gonadal hormones might represent a ‘sensitive window’ during which different sensory and motor systems in telencephalon and cerebellum can be seasonally re-shaped in both sexes.

1. MRI Data Processing

Whole-brain volume was manually delineated on the T2-weighted 3D anatomical RARE scan, which covered the entire brain, including telencephalon, diencephalon, mesencephalon and metencephalon. These volumes were used as a measure of brain size, used in further statistical analysis.

DW-images were prepared for voxel-based analysis using MRtrix3 version 3.0 (Tournier et al., 2012). We used an in-house algorithm to convert the Bruker 2dseq files to nifti files, which are compatible with other software programs such as SPM and MRtrix3. This step includes a signal scale correction, since the DWI data is acquired in three separate sequential scans. Furthermore, all DW-images were scaled by factor a 10 in 3 dimensions, to enable proper processing of small brains in software programs designed to process human data like SPM and MRtrix3. Voxel-based analysis requires that all images are spatially normalized to the same template, to enable voxel-wise comparisons. A simplified overview of the different MRI data processing steps is given in figure 1B-D.

1.1 DTI processing

Preprocessing the diffusion data for voxel-based analysis was performed using MRtrix3 (Tournier et al., 2012). First, diffusion gradient orientations were checked and automatically corrected to match the coordinate frame of MRtrix3 and ensure the best global ‘connectivity’ (Jeurissen et al., 2014a). Since the DTI data is acquired in three separate sequential scans, we further corrected for intensity differences between scans by rescaling the diffusion scans based on their b0 images. Preprocessing of the individual DW-images included the following steps: denoising (Veraart et al., 2016), correction for Gibbs ringing (Kellner et al., 2016), motion and distortion correction using FSL (Andersson and Sotiropoulos, 2016; Jenkinson et al., 2012), bias field correction using ANTS (Advanced Normalization Tool; (Avants et al., 2010)), creating an automated whole-brain mask for whole-brain extraction that were manually checked, upsampling to isotropic voxels of 1.75 mm. These preprocessed diffusion-weighted images were used to calculate individual diffusion maps (FA, MD, AD, RD). The transformation parameters derived from building the FOD template (see 5.5.2) were applied to the diffusion maps to warp them into the template space to perform voxel-based analysis. Next, these images were smoothed to double voxel size (3.5 x 3.5 x 3.5 mm³). Finally, all normalized diffusion maps were averaged to create an FA template that is used as a background to display the statistical results.

1.2 Fixel based analysis

For calculation of the fiber-based metrics, we followed the preprocessing steps as defined in (Raffelt et al., 2017). Fixel based analysis follows the same preprocessing steps as DTI processing up until the bias field correction. Apparent fiber density analysis differs from DTI analysis by the fact that it is related to the diffusion-weighted signal intensity within a given voxel. Therefore, global intensity normalization is performed to ensure robust voxel-wise comparison across subjects. Within this step, we used the default FA threshold of 0.4 to create an approximate WM mask, which is used to normalize the median white matter b=0 intensity across all subjects (Raffelt et al., 2012). For each image, a white matter response function was estimated for spherical deconvolution using the unsupervised Dhollander algorithm (Dhollander et al., 2019). Next, the average of all individual response functions was calculated and used for constrained spherical deconvolution to estimate FOD images (Jeurissen et al., 2014b). These FOD images were normalized to create an unbiased study-based FOD template, which involves linear and non-linear registration (Raffelt et al., 2011).

Next, the fixels in the FOD template are thresholded at 0.15, identifying the template white matter fixels to be included in further analysis. This threshold is lower than the default threshold of 0.25 for the human brain, as this threshold is too high for the songbird brain and excludes many of the genuine white matter fibers. We are aware that choosing a lower threshold comes with the risk of introducing noisy fixels, especially within grey matter, and take this into account in the interpretation of the results.

The estimated transformation parameters or warps of each subject to the template were used to transform the individual FOD maps into template space without FOD reorientation, so that the apparent FD can be estimated prior to reorienting the fixels. In the next step, we compare the fixels within a spatially matching voxel of the individual subject and the template white matter fixels, to identify which fixels correspond to each other and subsequently assign the corresponding FD value (Raffelt et al., 2017).

Next to FD, we also computed a fixel-based metric related to the macroscopic morphology in fiber bundle cross-section (FC). Fiber bundle cross-section (FC) information relies solely on the transformation parameters or Jacobian determinants generated during the construction of the population template, similar to other morphometry analyses like voxel-based morphometry (Ashburner and Friston, 2000). Morphological differences in the plane perpendicular to the fixel orientation could reflect differences in the number of axons, myelination, or the spacing between axons. For group statistical analysis, the FC values were logarithmically transformed to log FC to ensure that the data are centered around zero and normally distributed. Positive values indicate then expansion, whereas negative values reflect shrinkage of a fiber bundle relative to the template (Raffelt et al., 2017).

In each group, one scan failed due to excessive movement of the animal (subject ssw at SD12 and subject GoRR at LD4). In the male group, one animal died at the third time point (subject sblw). Of this animal we included MRI data of the first two time points.