Prairie vole reference brain and atlas generation for automatic c-Fos analysis.

(A) Generation of prairie vole (upper row) and mouse (lower row) reference brains using LSFM imaging. Top view of both prairie vole and mouse brains after perfusion. Cross-section of a coronal and horizontal views of prairie vole and mouse reference brains built after the co-registration of ∼200 brains per species. Prairie vole whole brain is 1.39 times bigger than the mouse (see Figure Supplement 1 and Supplemental File 1). 3D renderization of both prairie vole and mouse brains. Sagittal view of prairie vole and mouse whole brain fluorescent Nissl (NeuroTrace) staining imagined with LSFM. (B) Whole brain staining of both prairie vole (upper row) and mouse (lower row) brains registered to the reference brain. In red, boundaries of registered mouse reference atlas onto both prairie vole and mouse reference brains. NeuroTrace, TH and SST were registered and overlay to the atlases for validation. (C) Coronal sections of resulting prairie vole atlas after manual validation. (D) Overview of LSFM prairie vole c-Fos+ analysis pipeline. (left panel) Sagittal section of prairie vole c-Fos immunolabeling and imaged with LSFM. (center panel) Detail of two brain locations immunolabeled with c-Fos and overlay with resulting segmentation. (right panel) Area-based analysis of c-Fos+ cells. All c-Fos+ cells centroids are registered to the prairie vole reference brain and analyzed using the new prairie vole reference atlas. For each brain it is generated a voxelized representation of all c-Fos+ cells in the same prairie vole reference space and overlay with the reference atlas.

Study design and development of the prairie vole pair-bond.

(A) Schematic of the experiment design, where social behaviors (during 1h observations, black blocks) and IEG expression patterns are compared between mate pairs and siblings across time. (B) Continuous automated tracking of active social interactions (i.e., an index of social investigation and mating) for mate pairs (red lines) and sibling pairs (blue lines). (C) Time courses of specific social behaviors during 1h observation windows for IEG expression, with each row representing one pair. Social behaviors are overlaid onto self-grooming, and ultrasonic vocalizations (USVs, short red ticks) are overlaid onto all behaviors. (D) Plots showing group differences (mean ± sd) in individual activity level (velocity) and movement relative to the partner (net move, positive values indicate movement towards the partner). (E) Group differences (mean ± sd) in vocal behavior, proximity seeking, mating, and side-by-side contact. For (D) and (E) mate pairs are in red and sibling pairs in blue. Females are a lighter hue and males a darker hue for behaviors measured on an individual (rather than dyadic) level. T-tests were used to compare mates and siblings, and paired t-tests were used to compared female and male mates.

A brain-wide functional network is active during pair-bond formation.

(A) Map of bonding-associated voxels. Brightness indicates significance level of comparisons of hypothesized and null models (GLMs) to predict c-Fos cell counts. (B) Identification of significant and mutually exclusive ROIs, sorted by anatomical division and F-statistic. The significance levels of these models were determined with a permutation test, by comparing the observed F-statistic to a null distribution of F-statistics from shuffled data (n=10,000 permutations). Colored symbols for ROIs match their cluster group assignments in 3C and 3D. (C) Hierarchical clustering of chosen ROIs (n=68) and pairwise Pearson correlations of c-Fos cell counts. (D) Multi-dimensional scaling (MDS) coordinate space of the correlations between chosen ROIs. The most significant ROIs per cluster are labelled and their symbol size scaled by the F-statistic. Darkness and thickness of connecting lines reflect correlation coefficients. (E) Time course trajectories of total c-Fos cell counts within each cluster. Counts per cluster are scaled across samples and averages taken for each experiment group, with red lines for mates and blue for siblings (females-lighter, males-darker). Each cluster is given labels to summarize the most significant ROIs within them.

BST emerges as a central hub in the bonding network that is associated with mating success in both sexes.

(A) The first dimension of canonical correlation (CC) scores is compared across experiment group (mean ± sd). (B) Heatmaps represent correlation coefficients among CC scores, ROI c-Fos cell counts, and behavior measures. The full dataset is on top, and the bottom two correlograms are for female and male mates. For (B), (E) and (F), warm and cool colors represent positive and negative coefficients, respectively. Arrows mark the ROI (BSTpr) and behavior (ejaculation rate) with the strongest correlations to CC1. (C) BSTpr activation is compared across groups with a model comparison followed by a permutation test. (D) Successful mating events are shown across timepoints (mean ± sd). (E) Similarity (i.e., correlation) of bonding network activation is shown in female-male pairs. (F) Similarity in bonding network activation is shown in female-male pairs, using partial Pearson correlations to control for ejaculation rate. (G) Female-male pair similarity is shown for CC1 scores and BSTpr activation during bond formation. (H) Mating success is associated with BSTpr activation for both sexes during bond formation. Pearson correlations were used to assess these associations.

Working model for neural systems that shape stages of pair-bond development.

(A) Schematic of the network of regions identified in our study overlaid with regions known to be involved in rodent mating behavior (Pfaus & Heeb, 1997; Veening et al., 2014; Veening & Coolen, 2014). This network is proposed to be involved during the early stages of bond formation (e.g., 2.5 h timepoint). (B) Schematic of the network of regions identified in our study overlaid with regions known to be involved in prairie vole pair-bonding behavior (Walum & Young, 2018; K. A. Young et al., 2011). This network is proposed to be involved with the middle stages of bond formation when animals are engaged in prolonged mating and affiliative interactions (e.g., 6 h timepoint). (C) Schematic of the network of regions that are correlated between female and male mates at the 22 h timepoint in this study. These are identified from pairs of regions in which both sexes show high inter-individual similarity (Pearson correlation, r > 0.75). This network is proposed to be involved with the recognition and convergence of behavioral state in bonding partners. The network schematics in (A), (B) and (C) are adapted from the multi-dimensional scaling of correlations between region activity (Figure 3D). Light gray points represent regions not included in the proposed network. Line thickness and darkness between colored regions represent Pearson correlation r values between the connected regions.

Validation of the prairie vole reference atlas and brain-wide patterns of immediate early gene expression during pairing.

a. Mouse (left panel). Coronal sections composition of mouse reference brain overlay with atlas boundaries in red (left) and the same section of TH immunolabeling (right). Prairie vole (left column). Similar coronal brain sections of prairie vole reference brain and atlas boundaries in red (left) with TH immunolabeling (right). (center column) Same coronal section of prairie vole reference brain and c-Fos+ mean voxelization overlay with atlas boundaries in red. (left column) Same coronal sections of prairie vole NeuroTrace staining registered onto prairie vole reference brain. b. Coronal cross-sections (rostral to caudal) from female (top) and male (bottom) mating pairs are shown for the 2.5h timepoint group, with brightness corresponding to the average voxel c-Fos+ cell counts. c. A representative coronal slice that includes posterior BST is shown for mate pairs and siblings and separated by timepoint and sex. Brightness corresponds to average voxel c-Fos+ cell counts per group. d. Counts of c-Fos+ cells are shown for ROIs associated with pair-bond development, organized by hierarchal cluster. For each cluster, the total counts are shown on the right and the most significant ROI (highest F-statistic) is shown on the right. Counts are summarized by timepoint, partner type and sex (mean ± sd), and overlaid with counts from individual animals. Mate pairs are in red and sibling pairs in blue (female in lighter hues, males in darker hues). Test statistics are from an ANOVA to compare null and hypothesized general linear models. FDR q-values are from permutation tests to assess the likelihood of the observed test statistic as compared to a distribution of test statistics obtained from shuffled data. e. Representations of the first three dimensions of a multi-dimensional scaling (MDS) coordinate space based on Pearson correlations between c-Fos+ cell counts in brain regions (ROIs) associated with bonding. Each symbol represents an ROI and is colored based on cluster assignment. Darkness and thickness of connecting lines reflect correlation coefficients. f. This heatmap shows normalized connection densities (Knox et al., 2019) from projecting to ipsilateral target ROIs in the mouse brain. These ROIs are the same as, or larger divisions that contain, the ROIs found to be associated with pair-bond development in our analysis. g. The top histogram shows results of a permutation test to assess whether hierarchical clusters of chosen ROIs in our analysis mirror underlying anatomical connections. Observed connectivity reflects the overall average of cluster means, which is compared to averages from 10,000 iterations of shuffled data (i.e., target regions shuffled for each projecting region). The bottom histograms show the results of permutation tests for the mean connectivity within each cluster of ROIs.

Behavioral phenotypes during pairing and patterns of associate between brain-wide immediate early gene expression and behavior.

a. Long-term tracking of social states was informed by automated measures of the largest body area, pair activity, and pair distance. On the left, body area is plotted against dyad distance from an exemplar 24h mate pair video (10% randomly selected frames from 12 hours of white light). On the right, the density curve for body area reveals a basic threshold for when animals are separated or in physical contact. b. Traces of body area, video activity, and pair distance are plot for the first hour of interaction of the exemplar mate pair, along with automated assignments of behavioral states. c. Cumulative time spent in inactive social contact for up to 22h of cohabitation in mate pairs (red lines) and sibling pairs (blue lines). d. Group differences (mean ± sd) are shown for appetitive behaviors including nose-to-nose touching, anogenital investigation, and close follows. Group differences (mean ± sd) are shown for proximity and grooming behaviors including pair distance, allogrooming, and selfgrooming. Mate pairs are in red and sibling pairs in blue (female in lighter hues, males in darker hues). T-tests were used to compare mates and siblings, and paired t-tests were used to compare female and male mates. e. Hierarchical clustering of behavioral measures from Pearson correlations groups behavioral states and interactions into three main clusters involving close contact (e.g., huddling and allogrooming), mating (e.g., mounts, vocalizations), and appetitive behavior (e.g., approaches and follows). On the right, Pearson correlations are shown among behaviors in subsets based on partner type, sex, and timepoint. Warm and cool colors indicate positive and negative correlation coefficients, respectively. White indicates behaviors with no variation, meaning that coefficients were not computed (e.g., siblings did not exhibit mating behaviors). f. Schematic of canonical correlation analysis (CCA), where two sets of variables (x and y) are combined linearly to extract latent variates with the strongest correlations (CC1 stronger than CC2, and so on). This analysis outputs x and y scores for each latent variate for each animal subject in the dataset. In this dataset, x scores represent linear combinations of behavioral measures, and y scores represent linear combinations of ROI c-Fos counts. g. On the left, the relationship between x and y scores is shown for all study animals for the first three variates. On the right, the means of x and y scores for the first three variates are compared across partner type, sex, and timepoint (boxplots = median, 25%-75% quartile, and 95% confidence intervals). h. Heatmaps represent Pearson correlation coefficients, with brain regions on the y-axis (grouped into hierarchical clusters, see Figure 3 for ROI and cluster labels) and behavioral outcomes (during 1h observation windows) on the x-axis. Correlation heatmaps are split by partner type, sex and timepoint. Warm and cool colors indicate positive and negative correlation coefficients, respectively. White indicates behaviors with no variation, meaning that coefficients were not computed (e.g., siblings did not exhibit mating behaviors). i. On the top, heatmaps represent Pearson correlation coefficients between selected brain regions (ROIs) in mating pairs, with female data on the y-axis and male data on the x-axis. On the bottom, heatmaps represent Pearson correlation coefficients between selected ROIs in sibling dyads, where sibling A data are from animals placed in the left side of chamber during isolation, and sibling B data are from animals placed on the right side. ROIs are grouped by hierarchical clustering (see Figure 3). Warm and cool colors indicate positive and negative correlation coefficients, respectively.