Subjects were captive-bred sexually naïve male prairie voles descended from populations captured in southern Illinois. Voles were weaned on postnatal day 21 and housed with a same-sex conspecific in cages (29.2 L × 19.1 W × 12.7 H cm) containing corn cob bedding and crinkle paper nesting material with food (LabDiet Rabbit Diet 5321) and water ad libitum. Colony rooms were maintained on a 12L:12D photoperiod (lights on at 0600 h) and at a temperature range of 21 ± 1°C. Subjects were 70–100 days of age at the start of the experiment. All procedures were conducted in accordance with the National Institutes of Health Guide and Use of Laboratory Animals and the Institutional Animal Care and Use Committee at the University of Kansas.

The procedures used have been described elsewhere (Wang et al., 1997; Gobrogge et al., 2007). In brief, 40 total male subjects (10 per condition), aged 70–100 days were paired with sexually receptive females. Prior to establishing male-female pairs, females were estrogen primed with one daily subcutaneous injection of estradiol benzoate (EB) at 1 μg/100 μl sesame oil for 3 days consecutively, which successfully induced behavioral estrus. Pairs were housed together in plastic cages (29.2 L × 19.1 W × 12.7 H cm). Female partners were checked for vaginal plugs for the initial 4 days to ensure mating. On day 12, breeding pairs were transferred to larger cages (47.6 L × 26.0 W × 15.2 H cm) to establish residency. On the 15th day, male subjects were brought to the behavioral testing rooms and acclimated to their home cage for 5 min after their female partner was removed. Male subjects were randomly assigned to one of three social exposures for the resident social encounter test: (1) re-exposed to their familiar female, (2) exposure to an unfamiliar female, or (3) exposure to an unfamiliar male (both intruders age and size matched). Stimulus animals had a shaved patch on their back to distinguish them from the subjects. A fourth group of males was not exposed to any social stimulus to serve as the non-stimulated control. Behavioral interactions between the male subjects and the social stimulus were recorded for 10 min (between 0800–1100 h). Frequency and duration of affiliative and aggressive behaviors displayed by the resident was scored. This included olfactory investigation, side-by-side contact, huddling, allogrooming, attacks, bites, chases, defensive/offensive upright postures, offensive sniffs, threats, and retaliatory attacks. All behavior tests were manually scored by an experimenter blinded to experimental conditions using JWatcher (Blumstein and Daniel, 2007). The frequency and duration of affiliative and aggressive behaviors were analyzed by a one-way ANOVA (between subjects: condition). Pairwise comparisons between multiple groups were performed with Student-Newman-Keuls (SNK) post hoc analysis. Effects were considered significant at p < 0.05. After behavioral assessment, the social stimuli were immediately removed, and male subjects were housed alone for 2 h before being perfused (Winslow et al., 1993; Gobrogge et al., 2007). Female partners were sacrificed and dissected to determine pregnancy status. All females were determined to be 10–12 days pregnant based on previously validated methods in prairie voles (Aragona and Wang, 2009; Curtis, 2010).

Subjects were anesthetized with an i.p. injection of a ketamine / dexmedetomidine cocktail (75/1 mg/kg) then perfused through the ascending aorta with 15 ml of 0.9% saline, followed by 15 ml of 4% paraformaldehyde in 0.1 M phosphate-buffer (PB; pH 7.4). Brains were harvested, postfixed for 2 h in 4% paraformaldehyde at 4°C then stored at 4°C in 30% sucrose in PB for 3 additional days. Next, brains were cut into 30 μm coronal sections on a cryostat, and the slide-mounted tissue was stored in the −80°C freezer until processing for c-Fos immunohistochemistry.

A 1:4 series through each region of interest from each brain was removed from the −80°C freezer and allowed to come to room temperature for 5 min. Sections were rinsed in 1× Phosphate Buffered Saline (PBS, 0.1 M, pH 7.4) for 15 min (3 × 5-min) then with 1% NaBH₄ in 1× PBS at room temperature for 10 min. After rinsing in 1× PBS for 15 min, sections were blocked in 1× Powerblock (Universal Blocking Reagent 10×, Biogenex, Cat no. HK085-5K, Lot no. HK0850811) in 0.1% Triton X-100 for 30 min then incubated at room temperature for 30 min in rabbit anti-c-Fos polyclonal IgG (1:100, immunoStar, Cat no. 26209, Lot no. 1714001, RRID:AB_572267) in 0.05% TPBS (0.1 M PBS, 0.5% Triton X-100, pH 7.4) with 1% BSA then moved to 4°C for 48 h. Afterward, sections were rinsed in 0.05% TPBS for 5 min and then rinsed in 1× PBS for 10 min (2 × 5 min). Following, sections were incubated in 3% H₂O₂ in 1× PBS for 30 min then rinsed in 1× PBS for 15 min (3 × 5 min). The tissue was incubated at room temperature for 2 h in a goat anti-rabbit secondary antibody (1:150, Vector, Cat# BA-1000, Lot# ZF0809, RRID:AB_2313606). After, the tissue was rinsed in 1× PBS for 15 min (3 × 5 min) then incubated in VECTASTAIN ABC kit ((Peroxidase HRP), Cat# PK-4000, Lot# ZF1118) for 30 min. Next, the tissue was rinsed in 1× PBS for 15 min (3 × 5 min). Lastly, the tissue was stained with a DAB-Nickel Substrate Kit (Vector Lab, Cat# SK-4100, Lot# ZF0802) for 2.5 min, rinsed in 1× PBS for 15 min (3 × 5 min), and rinsed in ddH₂O for 1 min. Slides were dried and cover-slipped.

All tissue was imaged on a Leica DM6-B Microscope with a 10× objective. All regions were quantified bilaterally through ImageJ (Schindelin et al., 2012). Cell density (number of cells per mm³) of representative sections (3 brain sections) throughout each brain region was quantified by an experimenter blinded to experimental conditions and calculated. The cell counts for each condition was analyzed by a one-way ANOVA (between subjects: condition). Pairwise comparisons between multiple groups were performed with Student-Newman-Keuls (SNK) post hoc analysis. Effects were considered significant at p < 0.05.

A graph theory approach was used to analyze the functional network connectivity between all brain regions that were examined. Graph theory is a branch of mathematics and computer science that explores the patterns of connectivity between multiple vertices within a system. Prior research has employed graph theory to analyze large datasets across a wide array of disciplines including organic chemistry, public transportation, finance, etc., (Garcia-Domenech et al., 2008; Derrible and Kennedy, 2009; Quintero et al., 2013; Qiang and Fan, 2016). In addition, several studies have demonstrated functional connectivity associated with social recognition (Tanimizu et al., 2017), aggression (Teles et al., 2015), mating (Johnson et al., 2016; Kabelik et al., 2018), and positive social interactions (Rilling et al., 2018; López-Gutiérrez et al., 2021). We have used a graph theory approach to generate network models of functional connectivity between multiple brain regions. This approach allowed for an investigation of the most important nodes within the social decision-making network.

In most network models, nodes are grouped together into clusters based on shared characteristics. For example, a network model that reflects airport flight paths may group airports together into clusters based off of region. In this experiment, we employed a hierarchical clustering approach to unbiasedly cluster brain regions together based on shared patterns of coactivity. Interregional c-Fos correlations were used to calculate the Euclidean distances between each pair of brain regions for each behavioral condition. Distance matrixes were then arranged by both row and column according to their shared patterns of coactivity. Dendrograms were then generated to reflect the Euclidean distance between each brain region. The dendrograms were cut at half the height of each tree to cluster the brain regions into groups of nodes called modules. The dendrogram threshold of 0.5 was chosen because it has been used in prior c-Fos network modeling experiments (Kimbrough et al., 2018).

Both positive and negative edges were included in network analysis. Positive edges were included if two nodes had a Pearson’s R correlation >0.5. Similar functional network modeling experiments have used a range of edge thresholds spanning from R = 0.3 (Orsini et al., 2018) to R = 0.87 (Wheeler et al., 2013). We chose an edge threshold of R = 0.5 to ensure that all nodes within the network had at least one edge with one other node (Qi et al., 2014; Gao et al., 2019; Kimbrough et al., 2020). Negative edges were also included in cases where two nodes had a Pearson’s R correlation <−0.5. To our knowledge, this is the first experiment to measure c-Fos brain activity to included negative edges in network analysis.

To characterize individual nodes within the network, two metrics were used that reflect within-module connectivity and between-module connectivity. The within-module degree z score is a commonly used measure of the level of connectivity between one node and other nodes within the same module (Guimera and Nunes Amaral, 2005). Here, we developed a modified within-module degree variable that better reflects within-module connectivity for small network models. Our modified within-module degree z score (mWMDz) differs from the original in that it is calculated using the strength of the edges between multiple nodes instead of the number of edges between nodes. The mWMDz for node i, Eq. 1, is defined as:

r_{e_i} is the Pearson’s R values of edge e_i between node i and another node in the same module, s_i is the total number of nodes within module s, and K_i is the number of edges between node i and all other nodes within the same module. The mWMDz has a range of 0 to 1 with 0 reflecting no within-modular connectivity and 1 reflecting maximum within-modular connectivity. This range is achieved by using the absolute values of r_{e_i}, as our network analysis includes both positive and negative edges. If a node is in a module by itself, that node has a mWMDz of 0 by default K_i is the number of edges between node i and all other nodes within the same module.

The participation coefficient is a measure of how distributed the edges of a node are between multiple modules (Guimera and Nunes Amaral, 2005). Participation coefficient values range from 0 to 1 with 0 reflecting no participation and values approaching 1 reflecting greater participation. If a node only has edges with other nodes in the same module, it will have a participation coefficient value of 0. Alternatively, if a node’s edges are spread equally between many modules, that node will have a participation coefficient value close to 1. For participation coefficient, Eq. 2, K_is (between-module degree) is the number of edges (both positive and negative) between node i and nodes in module s. k_i (total degree) is the number of edges (both positive and negative) between node i and all other nodes within the network. The participation coefficient of each node, PC is then defined as:

All variables were calculated in python and then modeled in the Gephi 0.9.2 software package (Bastian et al., 2009).

Male prairie voles that were pair-bonded for 2 weeks exhibited different behavioral responses toward varying social stimuli presented. Male voles spent significantly more time interacting with the female stranger compared to the partner and stranger male (p < 0.05, Figure 1C). Further, males were predominantly affiliative when interacting with their partners and mainly aggressive when interacting with strangers (F_(2,27) = 9.98, p < 0.001; Figures 1B–D). The S-N-K test indicated that males exposed to an opposite-sex stranger displayed significantly more olfactory investigation, threat aggression (chases, offensive upright postures, lunges, etc.), and offensive aggression (attacks, bites, retaliatory attacks, etc.) compared to the partner. Male voles also displayed significantly more overall aggression (threat and offensive aggression) when exposed to a same-sex stranger compared to the partner (threat; F_(2,27) = 4.68, p < 0.02; offensive; F_(2,27) = 8.86, p < 0.001; Figures 1F–H). Our results of increased overall aggression toward the stranger conspecifics replicates previous results shown in Winslow et al. (1993) and Gobrogge et al. (2007). There were also group differences found in affiliative behaviors (F_(2,27) = 12.24, p < 0.0002, Figure 1E). Based on the S-N-K test, male voles exhibited significantly higher levels of affiliative behaviors toward their partner as compared to both the female and male stranger.

The immediate early gene, c-Fos is an established marker for neuronal activity and resulted in dense nuclear staining in eight of the thirteen brain regions associated with the SDMN (Chung and Brief, 2015). The resident exposure test induced c-Fos density increases in regional-specific and social stimulus-specific manners. First, group differences in c-Fos density were observed in the mPOA (F_(3,30) = 8.22, p < 0.0004; Figures 2M–P, 3D) and VMH (F_(3,30) = 4.45, p < 0.01; Figures 2Y–BB, 3G). c-Fos density differences were observed in the partner and same sex stranger in the VP (F_(3,30) = 5.16, p < 0.005; Figures 2A–D, 3A), LS (F_(3,30) = 5.19, p < 0.005; Figures 2E–H, 3B), and BNST (F_(3,30) = 5.64, p < 0.003; Figures 2I–L, 3C). Lastly, male voles exposed to same sex stranger induced a c-Fos density increase the AH (F_(3,30) = 5.52, p < 0.004; Figures 2Q–T, 3E), MeA (F_(3,30) = 3.04, p < 0.04; Figures 2U–X, 3F), and VTA (F_(3,30) = 9.91, p < 0.0001; Figures 2CC–FF, 3H). The S-N-K post hoc test indicated that male subjects exposed to one of the three social stimuli had increased c-Fos expression compared to baseline. Based on the S-N-K test, the increased c-Fos density is associated with the resident social exposure test. Our results show of c-Fos density increases in the BNST, mPOA, AH, MeA, and VTA shows similarity to previous research such as: Cushing et al. (2003), Curtis and Wang (2003), Gobrogge et al. (2007), and Wang et al. (1997).

To examine if various social exposures result in changes in brain activity and organization, we assessed the changes of neural coactivation of the SDMN and modular structuring as compared to the non-social exposed male voles. We examined the interregional coactivation of all regions using c-Fos correlations for each social exposure type. When first examined, there are clear differences between individual social and non-social exposure types (Figure 4). Specifically, the control voles showed higher resting-state coactivation throughout the network, where the voles interacting with varying social stimuli fewer positive correlations across the entire network, but instead, there were unique positive and negative coactivation patterns based on the social stimulus encountered (Figure 4A). When assessing coactivation patterns between social exposure types, the partner (Figure 4B) encounter increased positive coactivation patterns compared to both the opposite-sex and same-sex stranger (Figures 4C,D). Also, we see increased negative coactivation patterns for the same-sex stranger compared to the partner or opposite-sex stranger. Overall, these data suggest that distinct coactivation patterns are displayed in male voles depending on the social exposure and when re-exposed to their partner presents the strongest coactivation pattern between regions in the network.

Pearson’s correlation matrix allows for the observation of interregional coactivation as previously mentioned. However, this analysis only allows for the observation between two individuals regions. To examine the coactivation of a cluster of regions or regions that form modules, we used hierarchical clustering. Hierarchical clustering allowed for the observation and identification of the modular organization of the SDMN in non-social exposed control voles as well as the male voles exposed to one of three social exposures. Overall, the control voles displayed the least number of modules, indicating a general integration of connectivity in the SDMN during resting-state. By comparison, all social encounters remodeled the regional clustering and led to more restricted modules, with the exposure to the opposite-sex stranger generated the greatest number of modules including a number of regional isolates. We also found that there were no conserved modules across the social stimuli. Specifically, in control voles the network was organized into 3 distinct modules (Figure 5A). In the male voles re-exposed to their partner, we found that the network was organized into 5 distinct modules with one of the modules being an isolate (Figure 5B). Male voles exposed to an opposite-sex stranger demonstrated an increase, with a total of 8 modules and four of those modules are isolates (Figure 5C). Lastly, male voles exposed to a same-sex stranger displayed 5 total modules with a single isolate presented (Figure 5D). These data indicate that male voles present distinct modular structuring of the SDMN depending on the social exposure.

To further characterize our networks compared to the hypothesized SDMN, we used a graph theory approach. Refer to the methods section for specific calculations and terms. In short, we examined positive and negative connections, based on the Pearson’s correlation coefficient (threshold > ± 0.5R), between regions and modules displayed from the hierarchical clustering. We also determined the participation coefficient (strength of intermodular connectivity) and within-module degree Z-score (strength of intramodular connectivity) for all regions within the network.

The generation of these network models allows for better visualization of the functional connectivity of the SDMN during distinct social encounters beyond the correlational matrices and dendrograms (Figure 6; mWMDz and PC values are available in extended data Figure 6). Our network models suggest that male voles re-exposed to their partner had an increase in connectivity and organization as compared to the stranger conspecifics (Figure 6B). The partner network also demonstrates an overall increase in the participation coefficient and within-module degree Z-score as compared to the stranger conspecific networks. The opposite-sex stranger displays decreased connectivity and organization compared to the partner and same-sex stranger (Figure 6C). Both stranger networks demonstrated greater negative connections between modules compared to the resting-state or partner-evoked networks. This is most prominent in the same-sex stranger network (Figure 6D). Together, these data indicate that there are notable connective and modular differences within the functional connectivity networks based on the social stimulus presented.