Brain orchestration of pregnancy and maternal behavior in mice: A longitudinal morphometric study

Reproduction induces changes within the brain to prepare for gestation and motherhood. However, the dynamic of these central changes and their relationships with the development of maternal behavior remain poorly understood. Here, we describe a longitudinal morphometric neuroimaging study in female mice between pre-gestation and weaning, using new magnetic resonance imaging (MRI) resources comprising a high-resolution brain template, its associated tissue priors (60-µm isotropic resolution) and a corresponding mouse brain atlas (1320 regions of interest). Using these tools, we observed transient hypertrophies not only within key regions controlling gestation and maternal behavior (medial preoptic area, bed nucleus of the stria terminalis), but also in the amygdala, caudate nucleus and hippocampus. Additionally, unlike females exhibiting lower levels of maternal care, highly maternal females developed transient hypertrophies in somatosensory, entorhinal and retrosplenial cortices among other regions. Therefore, coordinated and transient brain modifications associated with maternal performance occurred during gestation and lactation.


Introduction
Motherhood is among the most transformative experiences in the lives of female mammals. While virgin females tend to avoid neonates, the end of the gestation period and the birth process lead to a behavioral switch characterized by an attraction towards infant cues, the expression of nurturing behavior and ultimately the establishment of infant bonding ( Lonstein et al., 2015 ;Numan and Insel, 2003 ). Decades of scientific research dedicated to the maternal brain have revealed a core neural circuitry that includes the medial preoptic area (mPOA) and the adjoining ventral part of the bed nucleus of the stria terminalis (BNSTv), which are highly critical for the onset of maternal behavior ( Bridges, 2015 ;Kohl et al., 2017 ;Numan and Insel, 2003 ). Func-Abbreviations: AC-PC, anterior commissure-posterior commissure; AMBMC, Australian mouse brain mapping consortium; AOB, accessory olfactory bulb; BNST, bed nucleus of the stria terminalis; CNS, central nervous system; CSF, cerebrospinal fluid; DARTEL, diffeomorphic anatomical registration using exponentiated lie algebra; df, degree of freedom; FA, flip angle; FLASE, fast large-angle spin-echo; FoV, field of view; MRI, magnetic resonance imaging; MOB, main olfactory bulb; mPOA, medial preoptic area; GM, gray matter; GMC, gray matter concentration; PVN, paraventricular nucleus of the hypothalamus; RARE, rapid acquisition with relaxation enhancement; ROC, receiver operating characteristic; ROI, region of interest; TE, echo time; TR, repetition time; VBM, voxel-based morphometry; WM, white matter. including structural plasticity through increasing neuronal soma size or astrocytic complexity within the mPOA ( Kinsley and Lambert, 2008 ), and changes in neurogenesis mainly in the main olfactory bulb (MOB)  but also in the mPOA/BNST in rodents ( Akbari et al., 2007 ). In human, regional morphological changes of gray matter (GM) within the parahippocampal gyrus, precuneus, cingulate, insula and frontal cortex have been observed in primiparous women using magnetic resonance imaging (MRI) ( Hoekzema et al., 2017 ). Additionally, olfactory cues coming from the neonate are integrated by the MOB and the accessory olfactory bulb (AOB) through an amygdalo-hypothalamic pathway, which is responsible for attraction/repulsion behavioral outcomes . Hence, the development of the maternal brain is dependent on both integration of external and internal cues acting through multiple brain pathways and regions to prepare the brain to gestation and motherhood. Nevertheless, the relationship between the brain rewiring over the gestation and lactation periods and the establishment of the maternal behavior is poorly documented.
To assess the dynamics of the maternal brain, a longitudinal MRI morphometric study over a complete reproductive experience was performed in mouse to investigate changes in the gray matter concentration (GMC) using voxel-based morphometry (VBM). VBM is a well-established and well-validated image analysis technique that provides an unbiased and comprehensive assessment of anatomical differences throughout the brain and has been successfully used to study GM changes within the mouse brain ( Ashburner and Friston, 2005 ;Keifer et al., 2015 ;Lein et al., 2007 ;Sawiak et al., 2014Sawiak et al., , 2013. However, available mouse brain MRI resources are often partial or provided in different spatial orientations or spatial resolutions ( Table 1 ). As an example, the Australian Mouse Brain Mapping Consortium (AMBMC), offers a high resolutive template and detailed atlases of the mouse brain including the cerebellum ( Ullmann et al., 2012 ), hippocampus ( Richards et al., 2011 ), diencephalon ( Watson et al., 2017 ) and cortices ( Ullmann et al., 2013 ). Unfortunately, segmentation of the MOB and the AOB and hindbrain is lacking, and this resource does not provide associated tissues probabilistic maps necessary for VBM. In other hand, the Allen Mouse Brain Common Coordinate Framework , is the most advanced mouse brain atlas ( Wang et al., 2020 ). This new atlas delimitates discrete structures within the thalamus, hindbrain, olfactory system and diencephalon and provides a full segmentation of cortical layers however, MRI template and brain tissues priors are still lacking for VBM investigations.
Thus, we combined these resources and associated probabilistic maps to emulate a complete resource dedicated to the mouse brain. Using these resources and a longitudinal VBM approach, we were able to assess the dynamic morphological changes of the brain during the whole reproductive period and demonstrate how these changes predict the quality of maternal behavior.

Animals
Twenty-three female RjOrl:SWISS virgin mice (8 weeks old; 20-25 g; Janvier Laboratory, Le Genest-Saint-Isle, France) were maintained on a 12-h light/dark cycle with access to food (standard chow) and water ad libitum . Animals were acclimatized 6 per cage to the housing facility for 7 days prior to manipulation. Females were randomly divided into two groups: a parous group ( n = 12), in which each female was exposed to a male (RjOrl:SWISS,8 weeks old; 20-25 g; Janvier Laboratory, Le Genest-Saint-Isle, France) for 5 days, became pregnant, and raised their offspring (litter size: 6 to 14 pups) until weaning at 21 days postpartum, and a control group ( n = 11), in which virgin females were not exposed to male mice. Each parous female was individually housed after male exposure. Control females were housed together in a separate room from parous females (two cages n = 5 and n = 6).
The MRI protocol was optimized to keep mice anesthetized for 2 h during each of the six acquisitions. During lactation MRI acquisitions, pups were kept under a heat lamp. One week after birth, maternal behavior was assessed as described in the behavioral section. All experiments were conducted in accordance with the local research ethics committee (APAFIS #6626-201002281145814V1) and are reported in accordance with the ARRIVE guidelines.

MRI acquisition
In vivo 3D MRI of the entire brain was performed three days before male exposure (baseline), at one week of gestation (early gestation), two days before the expected day of birth (late gestation), five days postpartum (early lactation), three weeks postpartum (late lactation) and two weeks after weaning (weaning). One female in the parous group and one female in the control group were scanned under similar conditions on the same day. Mice were anesthetized using isoflurane (2.5%; induction in O 2 /air mixture 1:1) (TEM-SEGA, F-33,600 Pessac, France) and then transferred and placed head first procubitus within an MRIcompatible cradle that incorporated a stereotaxic system dedicated for mouse head MRI, connected to a heater with circulating water to maintain body temperature and supplied with 1-2% isoflurane via a fitted mask. Respiration rate was recorded during all the experiments using an MRI-compatible monitoring system (PC-SAM model #1025; SA Instruments Inc., Stony Brook, NY, USA) and used to adjust the isoflurane rate to maintain a rate between 20 and 40 respirations per minute. After a recovery period of one hour, mouse returned to her pups. MRI studies were conducted at the Center de Biophysique Moléculaire d'Orléans and were performed on a 7T/160 mm PharmaScan spectrometer (Bruker Biospin, Wissembourg, France) equipped with an actively shielded B-GA09 gradient set, with 90-mm inner diameter and 300-mT/m gradient intensity. A 23-mm inner diameter Bruker birdcage coil with a cradle dedicated to a mouse head was used. Data acquisitions were performed on an Advance III console running ParaVision 5.1 software. T 2 -weighted images were acquired using a 3D fast large-angle spin-echo (FLASE) sequence which allows 3D brain mapping with a high resolution in a suitable time for in vivo acquisition ( DiIorio et al., 1995 ;Ma et al., 2008 ). Thus, the sequence with echo time (TE) = 20 ms, 1 repetition, acquisition matrix = 160 × 140 × 95, and a field of view (FoV) of 19.2 × 16.8 × 11.4 mm 3 , resulting in a final resolution of 120 μm isotropic voxels ( DiIorio et al., 1995 ;Ma et al., 2008 ). To obtain the FLASE sequence ( DiIorio et al., 1995 ), which is a specific sequence that is not included in the sequence package provided with ParaVision, the usual rapid acquisition with relaxation enhancement (RARE) spin-echo sequence was modified; in particular, the RARE-factor was fixed to 1 allowing a flip angle (FA) higher than 90°for the excitation pulse, while maintaining a 180°refocusing pulse ( DiIorio et al., 1995 ). Thus, T 2 -weighted images were obtained with a repetition time (TR) as short as 300 ms, 10 times lower than that needed for a classical T 2 -weighted spin-echo sequence. The sequence was optimized for acquisition in 1 h 28 min, with an isotropic resolution of 120 μm, a TR of 300 ms, an TE of 20 ms, an excitation pulse (FA) of 115°, with 1 repetition and with a matrix of 160 × 140 × 95 corresponding to a FoV of 19.2 × 16.8 × 11.4 mm 3 and contained the whole mouse brain.

Maternal behavioral test
One week after birth (P6), the maternal behavior of each female was assessed between 10:00 and 12:00 am using the pup retrieval test ( Keller et al., 2010 ). We chose to test mothers only once, in their home cages, to minimize the stress induced by the whole experiment which included six MRI sessions. Briefly, three pups were removed from the nest and placed at three different corners within the home cage. The latency to retrieve each pup and the time spent licking the pups, crouching in the nest over the pups and performing nonmaternal behaviors such as selfgrooming and digging were recorded over 15 min. Retrieval was defined Table 1 Comparison of mouse brain resources currently available in the literature (NA = not available). as the animal picking up a pup and transporting it to the nest. Crouching was defined as the animal assuming the nursing posture. Nursing and licking were permitted whether they took place in the nest. Digging was defined as a female out of the nest and digging the litter all around the nest. This version of the test permit to (i) evaluate the motivation of the dam through the measure of the retrieval of the first, second and third pups and (ii) evaluate the behavioural output of the dam after a moderate stressful stimulation. All videos were analyzed using BORIS software version 4.1.4.

K-means clustering
Clustering analysis of the behavioral data was performed using MAT-LAB Simulink 10b (The Mathworks, Inc., USA). Normality was verified and no outlier subjects were detected. To classify animals according to their maternal performance, a k-means clustering algorithm was used with crouching and digging times as behavioral markers ( Tables S1 and S2 ). Digging was chosen because it is indicative of high maternal stress ( Deacon, 2006 ;van den Brom et al., 2009 ). This algorithm iteratively grouped the animals by creating k initial centroids, assigning each animal to the closest centroid, iteratively re-calculating the centroids from the mean of its assigned animals and re-assigning the animals to each centroid until there were no more changes across iterations ( Magalhães et al., 2017 ). This clustering divided parous animals into a high maternal behavior group (with high crouching and low digging time) and low maternal behavior group (with low crouching and high digging time). We did not observe any significant difference in litter size between the two groups (low maternal, 9.5 pups ± 1.176; high maternal 8.5 pups ± 0.922; t -test p = 0.5185) and the use of all behavioral data did not modify clustering.

Mouse brain template and atlas building
For the MRI protocol, we developed a brain template and an atlas from the AMBMC brain template and the Allen Mouse Brain Common Coordinate Framework , respectively ( Fig. 1 ). First, we down-sampled the AMBMC template and its associated atlases and the Allen Mouse Brain Atlas and its associated Nissl images to a suitable resolution for MRI analysis (60-μm isotropic resolution; Fig. 1 , step 1 ). Then, all images were manually aligned to the anterior commissure/posterior commissure (AC/PC) axis, and the center of the images was defined relative to the AC ( Fig. 1 , step 2 ). The resulting template was then segmented into GM, WM and CSF probability maps using the unified segmentation approach ( Ashburner and Friston, 2000 ) of Statistical Parametric Mapping 8 ( SPM8 ) and the mouse brain priors provided by the SPMMouse toolbox ( Fig. 1 , step 3 ). In parallel, our T 2 -weighted anatomical images were realigned, coregistered, bias-corrected and normalized to our template. Using the SPMMouse toolbox, we also segmented the images as described above, and from these preprocessed images, we obtained a large set of 138 images for each tissue class ( Fig. 1 , steps 4-7 ). From these images, we built population-specific GM, WM and CSF priors. To build these priors, for each tissue class, we applied a diffeomorphic anatomical registration using an exponentiated lie algebra (DARTEL) approach, which is an automated, unbiased, and nonlinear template-building algorithm ( Ashburner, 2007 ) ( Fig. 1 , step 8 ). This new set of populationspecific tissue priors was used for both atlas building and final VBM preprocessing.
To normalize the Allen Mouse Brain Atlas to our brain template, we used the associated Nissl-stained images because (1) Nissl staining corresponds to the GM prior in terms of histology, and (2) this image was already coregistered to the atlas ( Fig. 1 , step 4 ). Therefore, we applied the segmentation function provided by SPM8 using the GM prior previously calculated from the Nissl image to generate the "Nissl2template " normalization matrix. We used this matrix to normalize the atlas to the template, while avoiding interpolation to maintain the label indices as integers ( Fig. 1 , step 7 ). Then, a visual inspection of each normalized label was carried out to assess whether the normalization process modified the position and volume of the structure too much. When necessary, holes were filled and labels were redrawn according to Paxinos and Franklin's atlas and using the FreeSurfer package. Finally, the olfactory bulbs and hind brain regions were completed, the corpus collosum and ventricles were drawn from the WM and CSF priors, and the cerebellum labels were replaced by the AMBMC cerebellum labels, which are more accurate. Finally, the atlas image was symmetrized (left-right). Our mouse brain template, priors and atlas were normalized within the same space and with the same final resolution (60-μm isotropic resolution), resulting in our final mouse brain atlas composed of a mosaic of 1320 ROIs covering the entire brain ( Supplemental Video 1 ).

VBM data preprocessing
Previously preprocessed normalized T 2 -weighted data were segmented into GM, WM, and CSF within SPM8 using the images of the population-specific priors ( Fig. 1 , step 9 ). Then, to produce a more accurate registration within each mouse as well as across all mice, a longitudinal VBM analysis was applied using the strategy described by Asami et al. (2012 ). First, a subject-specific template was created by the DARTEL algorithm using the previous tissue class images ( i.e. , GM, WM, and CSF maps) obtained from each mouse for the six time points. The DARTEL procedure releases individual-specific flow field maps, permitting the application of diffeomorphic normalization on images of each tissue class to spatially normalize each time point on a subject-specific template space. Each normalized tissue class image was modulated by the Jacobian determinant to account for the expansion and/or contraction of brain regions over time. Then, a population-specific template was created by the DARTEL algorithm using all subject-specific templates of the tissue class images. Here, the DARTEL procedure releases populationspecific flow field maps, permitting the application of diffeomorphic normalization of each animal onto the images of each tissue class. Finally, tissue class images were modulated by the Jacobian determinant, and the final modulated GM images were spatially smoothed with an isotropic Gaussian kernel with a 3-mm full-width at half-maximum and convolved with GMC images to create GMC maps ( Fig. 1 , steps 10-12 ).

VBM statistics and analysis
SPM8 was used to reveal the temporal and regional changes in the GMC maps. A second-level SPM analysis comprising a flexible factorial model, which is equivalent to a 2 × 2 mixed-model ANOVA with group as the between-subject factor and time point as the within-subject factor, was used to compare the control versus the parous groups and the low versu s high maternal behavior groups ( Barrière et al., 2019a ). The factors included in the analysis were subjects, group (control versus parous, or low maternal behavior versus high maternal behavior), and time points (baseline, early gestation, late gestation, early lactation, late gestation, and weaning). A brain mask was used to constrain the analysis within the brain and the significance threshold was set as p < 0.01 (t (126) = 2.356, control versus parous; t (60) = 2.39, low versus high maternal behavior, false discovery rate (FDR-corrected)), and the minimum cluster extent was set at 25 voxels. The results are presented on axial brain slice series generated by the Xjview SPM plugin. Corresponding surfacing results were produced with BrainNet viewer 1.6 ( Xia et al., 2013 ), allowing the generation of both brain meshes and brain plots to visualize data and create supplemental videos.

Postprocessing statistical analysis
To identify which regions each cluster correspond to, we used a combination of atlas-based and voxel-wise approaches. Clusters revealed by the flexible factorial analyses were binarized and the resulting masks were used as an inclusive ROI to extract GMC values from the corresponding regions of the GMC map using our mouse brain atlas and the Fig. 1. Processing of mouse brain templates and building an atlas from the AMBMC template and Allen Brain Atlas for data analysis and visualization. To create our resources, we used both AMBMC mouse brain template and the mouse Allen Brain Atlas and its associated Nissl images. ( 1 ) We down-sampled to a suitable resolution for MRI analysis (60-μm isotropic resolution) and ( 2 ) realign in the AC/PC axis. The resulting template was then segmented into GM, WM and CSF probability maps ( 3 ). These probability maps were used to segment all the images which have been previously normalized to the template ( 5,6,7 ). We obtained a large set of 138 images for each tissue class which have been used to build a population-specific GM, WM and CSF priors. using an exponentiated lie algebra (DARTEL) approach ( 8 ). This new set of population-specific tissue priors was used to segment again normalized T 2 images ( 9 ) for the final VBM preprocessing ( 10, 11, 12 ). To normalize the Allen Brain Atlas, we manually realign ( 13 ) and normalized the associated Nissl-stain mouse brain using the GM priors generated previously ( 14, 15 ). Both linear and nonlinear transformations have been applied to the Allen mouse brain atlas. Then, a visual inspection of each normalized label was carried out and, when necessary, redrawn according to Paxinos and Franklin's atlas. Finally, the olfactory bulb and hind brain regions were completed, the corpus callosum and ventricles were drawn from the WM and CSF priors, and the cerebellum labels were replaced by the AMBMC cerebellum labels. AMBMC = Australian Mouse Brain Mapping Consortium mouse brain template, AC/PC = anterior commissure/posterior commissure, CSF = cerebrospinal fluid, GM = gray matter, WM = white matter. REX plugin. As each cluster overlaps several regions in the atlas, the GMC data for each selected atlas region were extracted separately and then compiled and analyzed using GraphPad Prism 6.02 software to perform a region-based analysis. All regions with a size below 5 voxels were discarded. Regional GMCs were compared between groups and for each time point using a two-way ANOVA with repeated measures followed by a two-stage setup method of Benjamini, Kriegger and Yekutieli as recommended by the software (see supplementary tables). Maternal and nonmaternal behaviors (8 parameters) were compared using a multicomparison t -test with a false discovery rate (FDR) approach ( Q = 1%). Correlation analyses were restricted at two regions: entorhinal area ( Fig. 7 A ) and accessory olfactory bulb ( Fig. 7 B ). For these regions GMC data obtained at the late gestation time were correlated with both crouching and digging time using a parametric two-tailed Pearson test. Specificity and selectivity of the GMC data used for correlation to predict post-partum behavior were performed using the ROC curve method. Statistical significance was defined as p < 0.05 ( * ) for these analyses.

Mouse MRI atlas
VBM strategies require a template image and its associated priors of gray matter (GM), white matter (WM) and cerebrospinal fluid (CSF) for brain segmentation and normalization. In addition, a complete at-las of the mouse brain is mandatory for the identification of regions of interest (ROIs) highlighted by the VBM analysis. Given the limitations of available tools to thoroughly study GMC changes during the gestation and lactation periods in mice, we developed first a new set of resources using the AMBMC, an ultra-high-resolution template built from ex vivo brain images finely normalized within the same space ( Ullmann et al., 2013 ) and the Allen Mouse Brain Common Coordinate Framework ( Wang et al., 2020 ) ( Fig. 1 ). Our resources comprise the following four components: 1) a complete mouse brain template with a spatial resolution suitable for mouse brain analysis (60-μm isotropic resolution); 2) the corresponding GM, WM and CSF probabilistic maps for brain normalization together with a VBM analysis built from 138 T 2weighted images; 3) a complete mouse brain atlas derived from Paxinos and Franklin's mouse brain atlas ( "Paxinos and Franklin's the Mouse Brain in Stereotaxic Coordinates, Compact -5th Edition, " 2020 ) and composed of a mosaic of 1320 ROIs ( Fig. 2 A ); 4) a brain mesh permitting brain plot generation and data visualization ( Fig. 2 B and Supplemental Video 1 ). We visually inspected and carefully checked the results of the normalization process against the original coregistered atlas. Then, the labeled structures were reclassified and aggregated according to the brain regions to which they belonged (auditory, insular, temporal cortices, etc.), with respect to their anatomical topography (cortex, basal ganglia, etc.), tissue type (GM, WM and CSF) and hemisphere (left or right). Cortical structures were subdivided into functional ( e.g., primary and secondary motor cortices) or structural (agranular, Fig. 2. Details of the mouse brain template and atlas. ( A ) Coronal slices of the anatomical template of the mouse brain and the corresponding gray matter, white matter and cerebrospinal fluid probabilistic maps and the associated anatomical atlas (60-μm isotropic resolution). ( B ) Dorsal (left panel) and lateral (right panel) 3D representations of the anatomical mouse brain atlas. ( C ) Lateral views of the cortical areas after normalization of the Allen Mouse Brain Atlas to the AMBMC anatomical template. The cortex was segmented into cortical areas such as the cingulate, visual, auditory, frontal, visceral, insular, olfactory, retrosplenial, motor, somatosensory, temporal, entorhinal and gustatory systems. Each area was subdivided into secondary areas ( e.g., primary and secondary motor cortices) or structural areas ( i.e., agranular, dysgranular, agranular/dysgranular, granular and posterior agranular insular cortices). dysgranular, agranular/dysgranular, granular and posterior agranular insular cortices) areas ( Fig. 2 C ), as well as into different cortical layers ( Fig. 2 D ). Subcortical structures, as for example, the hypothalamus ( Fig. 2 E ) and the hippocampus ( Fig. 2 F ), were fully segmented according to Paxinos and Franklin's atlas (Paxinos and Franklin's the Mouse Brain in Stereotaxic Coordinates, Compact -5th Edition).

Morphometric changes occurred during the gestation and lactation periods
Next, we used our new resources to study the variations of GMC in mouse brain from the beginning of gestation until weaning. Using MRI T 2 -weighted anatomical acquisitions, we estimated the GMC maps which offer for each animal a global estimation of the GM. Longitudinal comparison of GMC between virgin mice (control group, n = 11) and mice who became pregnant and raised their young until weaning (parous group, n = 12) permits to highlight local modifications of GM during the whole reproductive cycle. A comparison of baseline and early gestation GMC maps between the control and parous groups did not reveal significant differences. However, at the end of the gestation period, significant increases in GMC were observed within several brain regions in the parous group compared to the control group ( Table S3 and Fig. 3 A ). Time course analysis revealed differences in GMC profiles between control and parous groups at the end of the gestation period, early in lactation and at the end of the lactation period. Specifically, GMCs within the mPOA and the BNST were consistently significantly higher in the parous group at that times. In addition, within the agranular insular cortex in the late gestation period and the early lactation period, GMC was significantly higher in parous group compare to control group ( Fig. 4 A ). During the early lactation period, we also found specific and significant increases in GMCs of the parous group within numerous brain regions ( Table S4 and Fig. 3 B ). Among these structures, the hippocampus ( dentate gyrus ; CA1, polymorph layer), amygdalar area and piriform area showed a transient increase in GMC at the early lactation time point that returned to baseline values at the end of the lactation period in parous group compared to control group ( Fig. 4 B ). In contrast, the caudate putamen, arcuate nucleus and paraventricular nucleus of the hypothalamus (PVN) showed significantly higher GMCs in the parous mice than in the control mice during the lactation period ( Table S5 and Figs. 3 C and 4 C ). Together, our data demonstrate that the late gestation period is associated with a pronounced increase in GMC in the mPOA/BNST, the core neural system of maternal motivation, lasting up to the late lactation period. Furthermore, the early lactation period is associated with increased GMC in other key maternal motivation areas in midbrain regions, including the hypothalamus, caudate putamen and amygdala. These GMC differences between both groups were no longer observed after weaning.

Morphometric changes during gestation predict the quality of maternal behavior
In the last part of this work, we evaluated whether these morphometric changes might reflect differences in maternal performance. Based on 15 min of behavioral observation during the pup retrieval test performed one week after birth, we evaluated the maternal performance of each mother by measuring the first, second and third pup retrieval times, pup-licking duration, crouching time, rearing time, digging time and self-grooming time ( Fig. 5 A ). Interestingly, we observed a large distribution of values for both crouching (284.2 s, SD ± 268.7 s) and digging (98.22 s, SD ± 159.3 s) durations within the parous group. Whereas crouching duration relates to maternal behavior, digging duration is widely recognized as a discriminative marker of stress-related behavior in rodents ( Deacon, 2006 ). Hence, we used crouching and digging durations to cluster animals using the k-means clustering procedure, thereby clustering parous animals into a high maternal behavior group (those with a high crouching time and low digging time) and a low maternal behavior group (those with low crouching time and high digging time).
A comparison of maternal performance parameters between the two clustered groups revealed significantly lower crouching and pup-licking times and a significantly higher digging time ( Fig. 5 B ) in the low maternal behavior group ( n = 6) than in the high maternal behavior group ( n = 6). Moreover, a comparison of GMC maps revealed both cortical and subcortical differences between the two clustered groups, mainly in the late gestation period ( Table S6 and Fig. 5 C ) but also during the early lactation period ( Table S7 and Fig. 5 D ). Indeed, transient increases in GMCs within the entorhinal area, lateral part of the orbital area, AOB, and medial preoptic area hat were observed at the end of pregnancy in the high maternal behavior group were absent in the low maternal behavior group ( Fig. 6 A ). In addition, the high maternal behavior group showed consistent higher GMCs in the hippocampus (field CA1, stratum oriens, field CA3), retrosplenial area and barrel field of the primary somatosensory cortex ( Fig. 6 B ) from the end of the gestation until the end of the lactation period.
Interestingly, using a receiver operating characteristic (ROC) analysis, we found that GMCs values within the entorhinal area and AOB at late gestation are reliable predictors for mouse maternal performance after birth. These GMC values significantly distinguished low maternal performance from high maternal performance postpartum (entorhinal area: sensitivity = 83.33, confidence interval (CI) = 44% to 99%; specificity = 83, CI = 44% to 99%; likelihood ratio = 5; Fig. 7 A ; AOB: sensitivity = 100, CI = 61% to 100%; specificity = 83, CI = 44% to 99%; likelihood ratio = 6, Fig. 7 B ). The GMC values of both the entorhinal area and AOB observed at late gestation were also significantly correlated with maternal behavior (crouching and digging times). These results reveal that the GMC differences in olfactory (AOB and entorhinal cortex) and mnesic (entorhinal area)-related brain regions occurring during the late gestation period significantly predicted the quality of maternal behavior.

Discussion
Using a new comprehensive neuroimaging resource dedicated to mouse brain, this longitudinal study reveals that pregnancy and lactation coincide with pronounced and transient cerebral changes. Transient increases in GMCs were observed in key regions controlling maternal behavior (mPOA, BNST, and PVN), as well as regions involved in emotions (amygdala), in motivation and reward (caudate nucleus, orbitofrontal cortex) and in mnesic functions (hippocampus). Interestingly, increase in GMC was also revealed in the insular cortex thought to link social and emotional skills. Moreover, we showed that females expressing high levels of maternal behavior had developed specific increases in GMCs in structures involved in olfactory (MOB and AOB) and somatosensory (somatosensory cortex) information processing, in memory (hippocampus, entorhinal cortex, retrosplenial cortex) and in reward and reinforcement (striatum). Interestingly, these hypertrophies were already significant at the end of the gestation period thus being predictive of the quality of maternal care ( Supplemental Video 2 ).

Implementation of new resources to support the analysis of mouse brain MRI data
The use of preclinical MRI is a target of growing interest for the study of brain structure and function in both healthy and pathological conditions. The use of advanced MRI techniques, coupled with the development of advanced animal models, is a powerful way to push new breakthroughs in the understanding of brain functioning and pathology. Herein, in the first step in our study, from the recent major advances in the development of brain mouse atlas, we generated a new set of neuroinformatic tools offering for the first time a complete resource dedicated to MRI studies of the mouse brain, namely, an accurate brain atlas (1320 ROIs), a high-resolution brain template and the associated  Orange and blue areas represent the gestation and lactation periods, respectively. Data are expressed as the mean ± standard error of the mean (SEM); two-way ANOVA followed by Holm-Sidak multiple comparisons test; * p < 0.05, * * p < 0.01 and * * * p < 0.001, compared with control mice.
GM, WM and CSF priors (60-μm isotropic resolution). The GM, WM and CSF probabilistic maps built and used in this study were calculated from 138 T 2 -weighted anatomical images, resulting in robust tissue class priors not only for VBM analysis but also for functional MRI and diffusion tensor imaging analysis in mice.
This comprehensive set of MRI compatible template and atlas for the mouse brain (Turone Mouse Brain Atlas and Template, TMBTA), allows a unified and standardized analysis of multimodal mouse brain MRI data and paves the way for the development of multicentric preclinical studies. Indeed, animal models deliver crucial information for the understanding of brain structure and function both in healthy and pathological conditions. Our template and mouse brain atlas were conceived to bridge the gap between basic and clinical neurosciences by providing to the preclinical neuroimaging community specific resources designed to be used in conjunction with the neuroinformatic tools and methodologies commonly used in human MRI studies. We anticipate that these resources will help neuroscientists to conduct their analyses of anatomical and functional datasets in a more standardized way, with the final goal of reaching more reproducible conclusions ( https://www.nitrc.org/projects/tmbta_2019 ).

Gestation and lactation periods induce strong but transient GMC hypertrophy
The establishment of accurate mouse MRI resources permits to study the variations of GMC longitudinally, in vivo and during the gesta-tion and lactation periods in female mice. We observed that several brain regions became transiently hypertrophic during pregnancy or in the lactation period until weaning. A set of structures comprising the core of the maternal circuit -especially the mPOA and BNST -displayed long-lasting hypertrophy that started at late gestation, culminated during the first week of lactation, and then disappeared at weaning. The initial GMC increase observed during the gestation period probably reflects changes induced by hormonal priming ( Brunton and Russell, 2008 ). Indeed, both the mPOA and BNST express a high number of steroid hormone and neuropeptide receptors ( McHenry et al., 2015 ;Tsuneoka et al., 2017 ). These factors are well known to trigger significant plasticity changes within the core maternal circuitry that are necessary for the preparation and adaptation of the brain to motherhood ( Bridges, 2015 ;Brunton and Russell, 2008 ). From parturition, the GMC differences observed during the whole lactation period highlight that mPOA/BNST receives a variety of sensory inputs from the pups, integrates that information with the females' endocrine status, and then projects to brain sites involved in socially-relevant motivation, affective state, and cognition ( Lonstein et al., 2015 ;Olazábal et al., 2013 ). Pup stimulation, electrolytic and neurotoxic lesions and local steroid hormone injections ( Bridges and Hays, 2005 ;Lee et al., 1999 ;Lonstein and De Vries, 2000 ;Miceli et al., 1983 ) have been shown to modify the intrinsic activity of these nuclei and consequently responsible for motivation and expression of maternal behavior. The mPOA is engaged throughout the postpartum period but differentially according to the needs of the developing pups. It has been shown that neurons Fig. 5. Distribution of animals according to the quality of their maternal behavior assessed with the pup retrieval test and brain morphometric. K-means clustering of parous animals to classify mice into low and high maternal behavior groups based on behavior during the pup retrieval test ( A ). Comparisons between the low and high maternal behavior groups revealed significant differences in crouching, pup-licking and digging times ( B ). Brain slices (left panel) and brain plots (right panel) comparing gray matter concentration (GMC) modifications and surface maps of GMC differences between females exhibiting low and those exhibiting high maternal behavior at the end of the gestation period ( C ) and early lactation period ( D ). Low and high maternal behavioral data were compared using a Student's t -test with post hoc corrections for multiple comparisons using an FDR approach ( Q = 1%) and are expressed as the mean ± SEM; * * p < 0.01 and * * * p < 0.001. SPM flexible factorial analysis revealed an interaction between low and high maternal behavior parous mice in the late gestation period ( A ) and early lactation period ( B ); voxel-level threshold p < 0.01, t (60) = 2.39, cluster threshold = 25 voxels. RS = retrosplenial cortex; vHippo = ventral hippocampus; S1 = primary somatosensory cortex; aOB = accessory olfactory bulb; OB = olfactory bulb; Pir = piriform cortex; Enth = entorhinal cortex; Str = striatum; OrbF = orbitofrontal cortex. Fig. 6. Longitudinal analysis in gray matter concentration (GMC) during the reproductive cycle in control females and females exhibiting low or high maternal behavior. Time-course comparison of GMC between the control (black dots and lines), low maternal behavior (pink dots and lines) and high maternal behavior (red dots and lines) groups revealed two types of time profile. GMC analysis in the entorhinal area, orbitofrontal area, the accessory olfactory bulb and the medial preoptic nucleus revealed an acute and specific increase in GMC values in the high maternal behavior group at the late gestation period ( A ). In contrast, GMC analysis in the ventral hippocampus, the retrosplenial area and the primary somatosensory area revealed an increase in GMC in the high maternal behavior group at the late gestation period, and this increase was maintained until weaning ( B ). Orange and blue areas represent the gestation and lactation periods, respectively. Data were compared using a two-way ANOVA followed by Holm-Sidak multiple comparisons test and expressed as the mean ± SEM; * p < 0.05, * * p < 0.01 when high maternal mice were compared with control mice and $$ p < 0.01, $$$ p < 0.001 when high maternal mice were compared with low maternal mice. Receiver operating characteristic (ROC) curves were estimated using a Wilson/Brown test with a 95% confidence interval. Correlations were estimated using a Pearson correlation test. ROC and correlation analyses were considered significant at p < 0.05. of the mPOA in late postpartum inhibits maternal responses allowing the changing expression and waning of maternal behavior across postpartum ( Pereira and Morrell, 2011 ). The sustained increase in GMC reported in late lactation could reflect the involvement of the mPOA to appropriately influence maternal behavior.
Our study highlights another set of brain structure that became hypertrophic only during the period from parturition to weaning -specifically, the PVN and arcuate nucleus of the hypothalamus. These changes illustrate the structural plasticity occurring in these regions. For instance, at parturition and during lactation oxytocin neurons of the PVN undergo dramatic neuronal, glial and synaptic changes such as an increase in size of the oxytocin neurons and an amplification of their synaptic input ( El Majdoubi et al., 2000 ). Oxytocin release at parturition facilitates the onset of maternal behavior by acting on the mPOA and is also important for maternal memory ( D'Cunha et al., 2011 ). Finally, lesions of the PVN disrupt the onset of maternal behavior ( Insel and Harbaugh, 1989 ;Numan and Corodimas, 1985 ;Olazábal et al., 2002 ). In the arcuate nucleus, dopaminergic cells are responsible for suckling induced prolactin release ( Freeman et al., 2000 ) and neurons projecting to the arcuate nucleus are involved in the maintenance of maternal motivation ( Cservenák et al., 2013 ).
Additionally, we report changes in GMCs found in olfactory related structures (main olfactory bulb, piriform cortex), somatosensory areas and auditory areas which reflect the multisensory control of maternal behavior. This finding is in accordance with a functional MRI study performed in rats which demonstrates that pup suckling is associated with increased neuronal activity within the midbrain, striatum and cortical sensory areas (somatosensory, olfactory and auditory cortices) ( Febo, 2011 ).
GMC variations within the hippocampus and entorhinal cortex, highlight the role of two essential structures involved in learning and memory processing during the reproductive period. Our data support evidence that the hippocampus undergoes profound neural changes during lactation. Indeed, lactating females have elevated spine densities in the hippocampus ( Kinsley et al., 2006 ) and show significant den-dritic remodeling in pyramidal neurons . Changes in hippocampal neurogenesis occurs during lactation and may support the enhancement of spatial memory necessary to foraging behavior in lactating females ( Duarte-Guterman et al., 2019 ;Eid et al., 2019 ;Kinsley et al., 2006 ;. Finally, our study also revealed a hypertrophy of the agranular insular cortex, which has never been reported in this context. The agranular insular cortex is a laminar part of the insular cortex and can be considered as a hub structure linking large-scale brain systems ( Gogolla, 2017 ). Indeed, the insula receives direct thalamic and somatosensory afferents carrying sensitive stimulations. In addition to its sensory afferents, the insula displays structural connectivity with the limbic system (basolateral, lateral and central amygdalar nuclei) as well as with the BNST, mediodorsal nucleus of the thalamus, lateral hypothalamus and perirhinal and lateral entorhinal cortices ( Gogolla, 2017 ). The insula also connects brain regions implicated in motivation and reward, such as the nucleus accumbens and caudate putamen ( Gogolla, 2017 ). Hence, our findings and the current literature suggest that before and after birth, the insular cortex may integrate and combine information from both external and internal stimulation and act as a relay between higher cortical and subcortical structures. Together, our results describe the dynamics of neurophysiological adaptation occurring in the brain from the early gestation period to weaning, thereby ensuring efficient maternal behavior and, by extension, the development of the offspring.

GMC modifications in the olfactory system at the end of the gestation predict the level of maternal behavior post-partum
All these transient modifications of GMC in parous animals indicate that various brain regions undergo considerable modifications from birth to weaning. Then, we sought to determine whether inter-individual variations in maternal behavior were associated with similar variations of GMC. Based on their behavioral performance during the pup retrieval test (crouching and digging times, Table S1 ), we used a k-means clustering strategy to divide maternal female mice into two groups displaying high versus low levels of maternal behavior. Even if we performed the pup retrieval test in the same conditions for each animal, the daily fluctuations of the maternal behavior could introduce a significant bias in the subsequent clustering. To validate our approach, we performed a negative clustering control using non-maternal behavior (rearing and grooming times, Table S2 ) demonstrating that our clustering procedure is specific to the maternal behavior and not due to the daily maternal behavior fluctuation. Using this strategy, we demonstrate that all mothers expressed similar maternal motivation since retrieving times were not significantly different. Nevertheless, the quality of maternal care was variable, half of the mothers spending higher time licking and crouching the pups and lesser time digging which indicates a highly adapted maternal behavior to a mildly stressful situation. The comparison between both groups revealed that several transient and long-lasting increases in GMCs were observed in the high maternal behavior group that were absent in the low maternal behavior group. Brain regions showing significant differences included the olfactory bulbs, somatosensory system, limbic system, especially the orbitofrontal area, and mnesic system, including the retrosplenial cortex, hippocampus and entorhinal area. Some of these structures are directly responsive to pup stimulation, and the observed dynamics (transient versus long lasting) may have been induced by different mechanisms during mother-offspring interactions (neurogenesis, gliosis, angiogenesis). For example, higher GMC values in the somatosensory cortex and olfactory bulbs in the high maternal behavior group potentially reflected increased suckling duration and proximity between the mother and pups, respectively. Regarding the differences in GMC in the ventral hippocampus, a large body of evidence reports that this region is critically involved in both anxiety and fear and it displays distinct properties of synaptic plasticity and intrinsic neuronal excitability from the dorsal hippocampus ( Bannerman et al., 2004 ). We hypothesized that the increase in GMC in the ventral hippocampus observed in high maternal dams reflects a local reorganization which could result in decreasing anxiety and fear. Such a general reduction in fear and anxiety does occur after birth for most females and has implications for many postpartum behaviors ( Lonstein et al., 2015 ). Mother rats genetically selected for their very low anxiety level are ineffective to respond to pups under challenging conditions ( Neumann et al., 2005 ). Thus, it could be that the differences in GMC in the ventral hippocampus between low and high maternal dams reflect differences in fear and anxiety resulting in differences in the quality of maternal behavior.
Strikingly, differences in GMCs were detected in the entorhinal area, orbitofrontal area, olfactory bulb, hippocampus, retrosplenial area and primary somatosensory area before parturition. These findings suggest that the maturation of these structures, probably through hormonedependent mechanisms, is a key determinant of the intensity of maternal behavior expressed during the lactation period. Using a ROC procedure, we found that GMC values of the entorhinal area and AOB at the end of gestation were significantly predictive of the maternal behavior postpartum. Interestingly, previous studies in mice reported an increase in cell proliferation during gestation in the subventricular zone, the neurogenic niche which provides newly generated neurons within the olfactory bulb (for review see ). These adult-born olfactory neurons are fully responsive to pup odor exposure ( Mak and Weiss, 2010 ) and are in part involved in some components of maternal behavior ( Larsen and Grattan, 2010 ;Sakamoto et al., 2011 ).
Hence, the correlations observed between GMC values in AOB at the late gestation period and maternal behavior performances may suggest that impairments of neural plasticity of this olfactory region would induce maladaptive neuroendocrine processing of the maternal brain at the end of the gestation period impacting maternal behavior performance. Taken together, our data provide the first potential imagingbased predictive biomarkers of the quality of maternal behavior and suggest the key role of the maturation of the olfactory system at the end of the pregnancy in the development of adaptative maternal behavior in mice.

Conclusion
Our study provides a new generation of neuroinformatic tools which will help basic neuroscientists to conduct structural and functional MRI investigations. Using these resources, we found that the development of the maternal brain is associated with substantial mesoscopic changes in critical regions. These modifications can be interpreted as cell size changes, neural or glial cell genesis/apoptosis, spine density or blood flow modifications ( Arentsen et al., 2017 ;Barrière et al., 2019b ;Kelly et al., 2015 ). As cellular and molecular modifications events are key for the adaptation of the brain to motherhood further molecular, cellular, and behavioral investigations must be performed to obtain a more precise view of the physiological mechanisms responsible for GMC variations. Nevertheless, the present study reveals a profound impact of motherhood onto the brain organization that results from two waves of changes, during the gestational and lactation periods (see Supplementary Video 2 where we display the temporal modification of brain). Importantly, we were able to distinguish the structural signature of subjects displaying distinct maternal performances paving the way for future functional investigations of these regions to improve our understanding of the maternal brain.

Declaration of Competing Interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Funding
The authors acknowledge the regional council of Center Val-de-Loire for funding this research through the IMACERVOREPRO grant (convention 201500104011 , 2015-2018) awarded to Matthieu Keller.

D.
A.B. contributed to the atlas and template building, data analysis, drafted and revised the manuscript. A.E. contributed to the imaging sequence troubleshooting, data analysis and contributed to the critical revisions. F.S. contributed to MRI sequence programming, imaging protocol setting, data acquisitions and contributed to the analysis and critical revisions. H.A., W.M., E.C., M.M., S.M. contributed to the study conception and design and contributed to the critical revisions of the manuscript. F.L. contributed to raise the funding, study conception and design and contributed to the critical revisions of the manuscript. M.K., is the principal investigator of the study, raised the funding, coordinated the project, revised and validated the manuscript.

Supplementary materials
Supplementary material associated with this article can be found, in the online version, at doi: 10.1016/j.neuroimage.2021.117776 .