Divergent neurogenomic responses shape social learning of both personality and mate preference

A major goal in evolutionary biology is to determine the mechanisms driving the origin and maintenance of reproductive isolation between species. Among closely related species, reproductive barriers are often reinforced by the development of behavioral suites that, in turn, are reinforced through cultural transmission and thereby sensitive to an individual's social environment (Verzijden et al., 2008, 2012a). Culturally transmitted homotypic mating preferences can rapidly promote assortative mating (Verzijden et al., 2005), whereas culturally transmitted personality traits, such as boldness, can accelerate local adaptation and reinforce ecological selection (Réale et al., 2010; Sih et al., 2012; Wolf and Weissing, 2012).

A growing body of work suggests that both learned and genetic variation in behavior is organized along a relatively small number of dimensions, collectively termed personality (Briley and Tucker-Drob, 2014; Réale et al., 2010). For example, individuals that are less likely to approach a predator are also less likely to approach a novel object – they are repeatably shyer across contexts, while other individuals are bolder. Personality, notably the shy–bold axis, often covaries with mating preferences within populations (Ingley and Johnson, 2014; Nosil, 2012; Schuett et al., 2010, 2011; Sommer-Trembo et al., 2016; Coleman et al., 2004; David and Cézilly, 2011; Both et al., 2005).

Covariance among behaviors is often the result of shared neural and hormonal mechanisms modulating suites of behavior across contexts, notably in the so-called social decision-making network (SDMN) (O'Connell and Hofmann, 2011, 2012; Simões et al., 2015). For example, variation in both personality and mate choice is associated with dopaminergic reward systems in the brain, and boldness and choosiness both covary with dopamine receptor expression in the SDMN (O'Connell and Hofmann, 2012, 2011; Cummings, 2015; Oliveira, 2013; Neckameyer, 1998; Thörnqvist et al., 2018). However, personality-related traits and mate choice decisions also rely on separate biological pathways. Boldness covaries with expression of stress-related genes in both dopaminergic (Thörnqvist et al., 2018) and serotonergic (Backström and Winberg, 2017) pathways of the brain, and immune response genes functionally link the immune system in the brain to anxiety-like behaviors and memory consolidation (Tonelli and Postolache, 2005; Marin and Kipnis, 2013). Mate choice measures, in contrast, are linked to a subset of pathways within the SDMN, including sex steroid hormone signaling (Davis and Leary, 2015; Sisneros et al., 2004), nonapeptide systems (Donaldson and Young, 2008) and synaptic-plasticity-linked genes such as neuroligins (Ramsey et al., 2014; Liang et al., 2015). Mating preference functions likely integrate these signaling pathways with social memory processing to make socially based decisions (Ramsey et al., 2014). Finally, mate choice depends on sensory perception and an individual's ability to integrate multimodal signals, and sensory gene expression can change markedly as a function of social and non-social experience (Ronald et al., 2012). By profiling behaviors and gene expression across social contexts, we can identify neurogenomic signatures pointing to brain mechanisms modulating suites of behaviors that are shaped in response to social ontogeny.

The goal of the present study was to apply an integrative approach to identify candidate molecular pathways implicated in the development of socially learned personality-related traits and mating preferences. To do so, we exposed female swordtail fish Xiphophorus birchmanni Lechner and Radda 1987 to one of three social environments (B-EXP: exposed to adult conspecifics; M-EXP: exposed to adults of the sister species X. malinche; NO-EXP: isolated from adult swordtails) throughout development and measured both conspecific mate preferences and boldness. We then conducted RNA-Seq on whole brains and olfactory epithelia to describe transcriptome-wide patterns in response to social exposure. To do so, we determined differences in the baseline regulation of: (1) our full transcriptomic data set, and (2) a subset of genes focusing on gene families related to the sensory periphery (odorant receptors) and hormone receptors, and other gene families implicated in social decision-making (O'Connell and Hofmann, 2012, 2011), personality (Thörnqvist et al., 2018; Backström and Winberg, 2017) and mate choice (Cummings et al., 2008; Cummings, 2015). We also used weighted gene co-expression network analysis (WGCNA) – a powerful approach for identifying clusters, or modules, of co-expressed genes and relating them to phenotypes of interest (Langfelder and Horvath, 2008) – to determine how these behaviourally relevant pathways covary with one another.

Swordtails rely on olfactory signaling for conspecific mate preference (McLennan and Ryan, 1999; Crapon de Caprona and Ryan, 1990; Fisher et al., 2006b; Verzijden and Rosenthal, 2011). Wild-caught female X. birchmanni normally prefer the chemical cues of conspecific males over males of their sister species, X. malinche. However, this preference is sensitive to both early (Verzijden and Rosenthal, 2011) and adult (Verzijden et al., 2012b) experience; female X. birchmanni learn to prefer familiar olfactory cues regardless of whether they are conspecific or heterospecific. By contrast, females in the sister species X. malinche show preferences for conspecific cues following exposure to X. birchmanni, but not to X. malinche (Verzijden et al., 2012b; Cui et al., 2017). These differences are mirrored by differences in boldness (Sih et al., 2004; Wilson et al., 1994; Toms et al., 2010). Xiphophorus birchmanni are bolder than X. malinche and show repeatable within-individual covariance among measures of boldness (Boulton et al., 2014; Johnson et al., 2015). In addition to sexual cues, swordtails attend to social information outside the context of mating (Wong and Rosenthal, 2005; Coleman and Rosenthal, 2006), suggesting that boldness might also be sensitive to social experience. In the present study, we show that social experience simultaneously drives intraspecific differences in both learned personality and mating preferences, and that the neural mechanisms driving the development of these behaviors are themselves dependent on social upbringing.

We pooled offspring of 15 X. birchmanni females (collected from Río Coacuilco at Coacuilco, Hidalgo, Mexico) within 3 days of birth and randomly assigned them to one of three treatment groups (30 juveniles/tank, three replicates): (1) B-EXP: visually and chemically exposed to two males and two females of adult X. birchmanni from the Río Garces locality; (2) M-EXP: exposed to two males and two females of adult X. malinche from the Chicayotla locality (Culumber et al., 2011); and (3) NO-EXP: controls that did not receive adult stimulus exposure. Exposure treatments were performed at 23°C, 12 h:12 h light:dark cycle in adjacent, though visually isolated, 120-liter aquaria in which adults and juveniles were divided by a transparent, perforated Plexiglas board, which allowed for transmission of both visual and olfactory cues (Verzijden and Rosenthal, 2011). Sufficient shelter was provided to both adults and juveniles, and water was continually refreshed through a flow-through system. Juvenile males were removed from aquaria upon the first sign of maturation (thickening of the anal fin to form the gonopodium).

Once all females had matured (approximately 9 months), we tested female preference for olfactory cues of X. birchmanni and X. malinche. In addition, we recorded measures of shy–bold behavior commonly used in swordtails and other fish (Boulton et al., 2014; Dingemanse et al., 2007). After all behavioral trials, we returned females to their respective treatment for an additional 2 months prior to tissue collection. We chose this time frame to minimize short-term effects from mating preference trials, which have previously been studied in other systems and contexts (Bloch et al., 2018; Wong et al., 2012; Wong and Cummings, 2014), and focus our analysis instead on the effects of social exposure on baseline gene expression profiles.

We obtained approval for experiments from Texas A&M University's Institutional Animal Care & Use Committee (AUP 2015-0112), and all experiments complied with current state, federal and local laws in the United States and Mexico.

We tested female preference for conspecific versus heterospecific male odors following a well-established protocol (McLennan and Ryan, 1999; Verzijden and Rosenthal, 2011; Cui et al., 2017). Briefly, a peristaltic pump dripped olfactory cues of model X. birchmanni and X. malinche males on either side of a trial tank. Briefly, to produce the olfactory cues, 20 liter aquaria were thoroughly rinsed with Alconox and a 1:1 mixture of hydrogen peroxide, followed by rinsing with carbon-filtered water four times. Groups of four males of X. birchmanni and four X. malinche were separately placed in 18 liters of carbon-filtered water and visually exposed to five females from their own population in adjacent tanks for 4 h. None of the males or females used to produce olfactory cues were used as exposure models. All preference trials were conducted with the same stimulus water. The preference trial tanks (75×19×20 cm filled to a depth of 15 cm) were opaque on all sides. A small shelter was provided in the middle of each tank. Female position during the trial was recorded using an overhead camera connected to Viewer (Biobserve GmbH, Bonn, Germany) recording software. The tank was equally divided along its length into three virtual zones: two preference zones on either end of the tank, with the middle defined as a neutral zone. Twenty minutes before each trial, the focal female was introduced to the testing tank for acclimation. Then, stimulus water started dripping on both ends of the tank via computer-controlled peristaltic pumps (VWR) until the end of the trial at a rate of approximately 5 ml min⁻¹. Upon initiation of the cue pumps, we allowed 5 min for the focal female to visit both preference zones. If the subject failed to do so, she was defined as unresponsive. Starting at the moment the subject visited both preference zones, the time in each zone was recorded for a total of 5 min. To account for potential side biases, we tested each female twice in succession, with the location of each cue switched. We summed the association times from both trials for data analysis. At the time of video scoring, all those who scored were blind to focal individals' treatment and replicate assignments. We tested a total of 30 B-EXP (22 responsive), 30 M-EXP (21 responsive) and 31 control females (27 responsive).

For association time datasets, we used one-tailed Wilcoxon signed-rank tests to detect differences in mean association time between the two stimuli for each group, with the a priori hypothesis that B-EXP and M-EXP female X. birchmanni would show preferences for olfactory cues of their exposure models as in Verzijden and Rosenthal (2011). We conducted a random-effects ANOVA on net preference (time with X. birchmanni – time with X. malinche in seconds) using the nlme package (https://cran.r-project.org/web/packages/nlme/index.html) in R version 3.6.0 (https://www.r-project.org/) to test whether exposure experience significantly affected mate preference. Replicate was included as a random effect, and exposure treatment was the fixed effect. Unpaired t-tests assuming equal variance were then used to test for differences in net preference between groups.

In this study, we extracted common measures of shy–bold behaviors from olfactory preference trials. Specifically, we measured total distance traveled (D_t) as a measure of exploratory behavior. Briefly, every 20 frames (2 s) we measured an individual's location in x–y space (origin placed at the bottom left of the arena). We then calculated D_t as the sum of all linear distances between successive data points. In addition, we measured the amount of time females spent within the shelter provided during preference trials (T_s). We chose these measures as they have each previously been suggested to be reliable indicators of boldness (Boulton et al., 2014; Dingemanse et al., 2007). We used unpaired t-tests to test for differences between exposure groups in both D_t and T_s. To examine the relationships among personality-related traits, we conducted a random-effects ANOVA on T_s of all females using the nlme package, including exposure and D_t as fixed effects and replicate as a random effect. Separate random-effects ANOVAs were conducted for D_t and T_s to examine the relationships between each personality-related trait and net preference of responsive females within exposure groups. Exposure and net preference were treated as fixed effects and replicate as a random effect. All analyses were conducted in R.

We randomly selected five females from each replicate of each treatment group (45 total samples). Females were retrieved from treatment tanks and immediately euthanized by placing them in a 250 ml beaker with an overdose of MS-222, and then decapitated, as previously done in X. malinche (Cui et al., 2017). In total, this process took less than 2 min per individual. For RNA-Seq analyses, we preserved whole heads in RNALater solution, placed at 4°C overnight, then dissected whole brains, optic nerves and both nares, including olfactory epithelia, and stored them at −80°C until RNA extraction.

We extracted total RNA from the above tissue using a standard Trizol reagent (Life Technologies) protocol following the manufacturer's instructions, and quantified and assessed the RNA for quality on a Bioanalyzer 2100 (Agilent Technologies). Only those samples with RNA integrity numbers >7 were used for further analyses (B-EXP, n=11; NO-EXP, n=13; M-EXP, n=12). We constructed sequencing libraries with Illumina's TruSeq RNA Library Preparation Kit v2 (Set A), following the manufacturer's protocol. Pooled libraries were sequenced on three Illumina HiSeq 2500 lanes (125×125 paired-end reads). We used cutadapt (Martin, 2011) to clean the raw reads (leading, trailing and sliding window quality ≥20 PHRED scale). Only reads >30 bp after filtering were analyzed. We then mapped trimmed reads to the previously described X. maculatus reference genome (Schartl et al., 2013) using Tophat 2.0 (Ghosh and Chan, 2016). First, we mapped pooled reads from all individuals to obtain a comprehensive alternative junction list. We then mapped reads for each individual sample separately guided by this junction list. We discarded reads with gap lengths greater than 1, or with more than 5 mismatches or edit distances.

We downloaded gene models (v. 82) for X. maculatus from Ensembl (Yates et al., 2015). We counted the number of reads mapping to each gene using the python package htseq-count (Anders et al., 2015; strand specific: no; mode: union; counted feature: exon), requiring a mapping quality of 20. We imported these raw counts into the DESeq2 package in R (Love et al., 2014) for differential expression analysis. We visualized gene expression profiles of individuals by conducting a constrained correspondence analysis on DESeq2 normalized gene counts using the vegan package in R (Oksanen et al., 2007). We performed a transcriptome-wide analysis, following the DESeq2 default parameters. We excluded genes that did not surpass a 0.5 counts per million threshold in at least four individuals from analyses. To control for differences among tanks, replicate was included in the design formula. We defined significance as genes differentially expressed between treatment groups at P<0.05 after Benjamini–Hochberg adjustments.

To compare how the social environment affects transcriptional regulation in female X. birchmanni and the sister species X. malinche, we identified shared genes that were differentially expressed between conspecific- and heterospecific-exposed females in both species. To identify the functional relevance of these genes, we used the human Allen Brain Atlas (Dong, 2008) to identify where each gene within our list is constitutively expressed within the brain. We conducted Fisher's exact tests in Enrichr (Chen et al., 2013) to identify brain regions that were enriched within our overlapping gene list. We then conducted literature reviews within the Web of Science Core Collection to identify each gene's functional role within the brain. Relevant results, along with a list of citations, are appended to Table S3.

In addition to a general, transcriptome-wide analysis, we created a list of a subset of gene families previously shown to be implicated in SDMN (X. maculatus annotated genes related to vasopressin, oxytocin, androgen, progesterone and dopaminergic signalling), mate choice decision-making (neuroligins and tyrosine hydroxylase), olfaction (odorant and vomeronasal receptors) or stress-coping (genes related to neuropeptide Y or dopaminergic signalling). We assessed differential expression between B-EXP and M-EXP individuals with a P<0.05 Benjamini–Hochberg adjustment, and noted WGCNA module assignments for each gene to determine how the genes from the different molecular pathways covary with one another. To visualize the differentiation between exposure treatments in SDMN/stress-related genes and olfaction-related genes, we conducted NMDS through the vegan package (Oksanen et al., 2007) in R, using normalized transcript counts as input. PERMANOVA was conducted (1000 permutations) to determine whether exposure had a significant effect on the expression profiles of these gene subsets: (1) SDMN/stress-related genes and (2) odorant/vomeronasal receptors (see Table S7 for full list of genes). We further examined the six highest influencing genes (based on highest absolute values of loading scores) for gene expression differences between exposure groups.

Systems genetics uses the connectivity of genes to describe the relationship between the transcriptome and a trait of interest. To identify modules of interconnected genes that correlate with social treatments, we used WGCNA (Langfelder and Horvath, 2008) on transformed count data for genes that passed an initial 0.5 CPM threshold filter. We excluded two outliers (one B-EXP and one M-EXP sample) from subsequent WGCNA analyses, based on a standardized connectivity threshold of 2.0 standard deviations from mean connectivity. We then reanalyzed sample network connectivity among the remaining samples. For all possible pairs of variable genes, we calculated Pearson correlation coefficients across all samples. We created an unsigned matrix and adjusted the soft-threshold value to ensure a scale-free topology (β=12), thereby creating a weighted network. Within this topological overlapping network (Yip and Horvath, 2007), genes were hierarchically clustered, and modules were identified based on the degree of similarity among genes (Langfelder and Horvath, 2008). We used a merging threshold of 0.2 (mergeCutHeight=0.2), with a minimum module size of 30 (minModuleSize=30) and a mean connectivity threshold of greater than or equal to 0.7 (minkMEtoStay=0.7). We used default parameters for the rest of the analyses. We then correlated module eigengene values for a given module across samples via Pearson's correlation, and identified modules differentially expressed between groups at P<0.05.

To visualize WGCNA genetic covariance results among modules significantly associated with social experience, we exported final co-expression networks to Cytoscape (Shannon et al., 2003). We attached information on module membership, whether a gene was differentially expressed, and whether a gene was a member of a gene family related to the SDMN pathway or stress-coping to the network as metadata so this information could be visualized. To identify SDMN/stress-related genes that may be considered to have more central functions within and across networks, we used the Network Analyzer tool within Cytoscape and ranked genes based on three scores: (1) betweenness centrality, which reflects the amount of control a given gene has on the interactions of other genes within a selected network (Yoon et al., 2006), (2) degree, which is the number of connections a given gene has with all other genes within a network, with more central genes having higher degree scores than peripheral genes, and (3) clustering coefficient, a measure of how connected a given gene is to the rest of the selected network of genes (Shannon et al., 2003). We conducted this analysis examining: (1) the 500 most central genes among the SDMN/stress, olfaction and vision modules, and (2) all genes directly connected to annotated olfactory receptor genes. We conducted Fisher's exact tests to determine whether gene list 1 or the 50 most central genes in gene list 2 were significantly enriched for annotated SDMN/stress-related genes (see Table S4 for list of genes).

To determine whether particular functional categories might be implicated in developing socially learned behaviors, we performed gene ontology (GO) enrichment analysis on: (1) the three lists of differentially expressed genes between pairwise comparisons of the three exposure groups (B-EXP/M-EXP, B-EXP/control and M-EXP/control), and (2) WGCNA modules that were differentially regulated (P<0.05) between exposure treatments (B-EXP/M-EXP). We used the annotated X. maculatus genome to assign Human Genome Organization gene symbols to each gene. We included all genes that passed coverage filtering in DESeq2 as part of the gene universe, and we tested for significant enrichment (FDR <0.05) of different biological processes by comparing the gene universe with gene lists using the PANTHER Classification System (release 20160321) (Mi et al., 2013). Briefly, genes are organized into families and subfamilies according to sequence homology and functional similarity. They are then assigned GO terms (‘GO biological processes complete’). Gene lists are compared with the gene universe to find GO terms that are statistically overrepresented or underrepresented using a binomial test. We used Revigo (Supek et al., 2011) to visualize GO categories clustered by semantic similarities (SimRel; Schlicker et al., 2006).

Early social experience had a significant effect on olfactory preference (F_2,65=5.42, P=6.7×10⁻³; Fig. 1A), with B-EXP and M-EXP females preferring the odors of males they experienced (Wilcoxon signed-rank tests, B-EXP, Z=2.65, P=4.02×10⁻³; M-EXP, Z=−1.65, P=0.049). NO-EXP females failed to show preferences for either cue (Z=0.120, P=0.452). B-EXP females had significantly greater preferences for X. birchmanni cues than both M-EXP (P=2.0×10⁻³) and NO-EXP females (P=0.034).

Exposure type had a considerable role in shaping female boldness, as heterospecific-exposed female X. birchmanni spent significantly more time within the provided shelter than conspecific-exposed females (M-EXP: 41.6±11.7 s, B-EXP: 15.9±4.2 s, P=0.045; Fig. 1C). No-exposure females were intermediate (NO-EXP: 17.6±5.3 s, P>0.05 from both M-EXP and B-EXP). Furthermore, heterospecific-exposed female X. birchmanni traveled significantly shorter distances (M-EXP: 456.6±50.9 cm) than both conspecific-exposed (B-EXP: 685.3±94.2 cm, P=0.038) and isolated females (NO-EXP: 701.2±71.8 cm, P=7.53×10⁻³; Fig. 1B). Personality-related traits did not significantly covary with one another (F_1,84=1.14, P=0.29) nor with female mating preference in our full models (D_t: F_1,61=1.67, P=0.20; T_s: F_1,61=1.00, P=0.32).

Social exposure had a large effect on gene regulation. After coverage filtering with DESeq2, 19,126 genes were retained, of which 17,715 had annotated gene symbols. Transcriptomic results mirrored behavioral results, with B-EXP and M-EXP females showing the greatest level of transcriptomic separation. NO-EXP females showed high variance in gene expression profiles, resulting in relatively little separation from either exposure treatment (Fig. 2A). After false discovery correction at FDR=0.05, 919 genes were differentially expressed between B-EXP and M-EXP females (Table S1). GO enrichment analysis revealed many of the genes upregulated in B-EXP females to be related to immune response – namely, the activation and regulation of T cells and leukocytes and tumor necrosis factor production, as well as behavioral response to stimuli and ribosome assembly. M-EXP females showed upregulation in genes related to neurogenesis and synaptic plasticity (Table S2, Fig. 3A,B). Sixteen of these genes were also differentially expressed between B-EXP and M-EXP X. malinche females in a similar study (Cui et al., 2017). Literature review and analysis of this gene list with the Allen Brain Atlas revealed that the majority of these genes (12 of 16) are associated with functional roles such as neuronal plasticity, learning and addiction (Dong, 2008). This overlapping list was enriched for genes associated with the lateral septal nucleus (Fisher's exact test, adjusted P=0.043), a region characterized as a functional connection between the social behavior network and the mesolimbic reward system with known roles in social recognition and the evaluation of stimulus novelty (O'Connell and Hofmann, 2011) (Table S3). By contrast, there was only a very weak signal of differential expression between NO-EXP and either B-EXP (six genes with significant differential expression) or M-EXP (two genes).

In addition to examining transcriptome-wide data, we assessed differential expression between B-EXP and M-EXP females in key gene families previously shown to be implicated in social decision-making (O'Connell and Hofmann, 2012), mate choice (Cummings, 2015; Wong and Cummings, 2014; Wong et al., 2012; Bloch et al., 2018) and stress-related personality traits (Backström and Winberg, 2017; Thörnqvist et al., 2018). Interestingly, we found that M-EXP females were characterized by greater expression of several genes related to the signaling pathways of dopamine, serotonin, oxytocin and progesterone, as well as neuropeptide receptor-related GPCRs, while B-EXP females mostly showed upregulation of olfactory receptors and MHC-I-like genes (Fig. 4, Table S4).

WGCNA revealed social experience-associated gene modules only between B-EXP and M-EXP treatments. Comparing B-EXP and M-EXP gene expression profiles, genes were grouped into 12 modules of similarly coregulated genes (Fig. 2B,C), three of which were differentially regulated between B-EXP and M-EXP females (Fig. 2D). We refer to these modules according to their functional roles [SDMN/stress (blue) module, 2768 genes, P=0.02; vision (purple) module, 86 genes, P=0.02; olfaction (yellow) module, 952 genes P=0.01; Table S5]. GO enrichment analysis revealed that the SDMN/stress module, upregulated in M-EXP females, is largely composed of genes with GO terms similar to those found in the differential expression analysis (Table S2). Furthermore, a large majority of genes related to social and mate choice decision-making are found in this module. The vision module, upregulated in B-EXP females, consists of genes with roles in visual perception (Fig. 3C, Table S2). The olfaction module, upregulated in B-EXP females, consists of genes with roles in ribosomal activity as well as a wide variety of GO terms, including the majority (54%) of annotated swordtail odorant/vomeronasal receptors (Fig. 3D, Table S2). Furthermore, all show the same trend of relative downregulation in M-EXP females relative to B-EXP and NO-EXP females (Fig. 4B,D, Table S4).

Based on the among-gene covariances calculated in WGCNA, we examined the connectivity within and among significantly differentially regulated modules between B-EXP and M-EXP females (SDMN/stress, vision and olfaction modules) to determine how these modules are neurogenomically linked. These modules were significantly enriched for SDMN or stress-related genes within our annotated gene list (Fisher's exact test: P=0.0217; Table S6). Four of these genes, all members of the SDMN/stress module, were among the 500 most central, or ‘hub’, genes across all three modules, with the neuroligin nlgn2b and neuropeptide Y receptor npy8ar both showing considerably high degrees of connectivity [in order from highest to lowest ‘betweenness centrality’ scores (Shannon et al., 2003; Assenov et al., 2007): nlgn2b, npy8ar, npffr2 and gabrg2; Fig. 5, Table S6].

As chemical communication is the primary driver of species-specific mate choice within the X. birchmanni–X. malinche system, we also examined the network connections within odorant and vomeronasal receptor genes (Fig. 6, see Table S4 for list of odorant/vomeronasal receptor genes). Similar to our analysis across all three differentially expressed modules, the 50 most central genes within this gene list were enriched for SDMN/stress-related genes (Fisher's exact test: P=3.20×10⁻⁴). The GABA-A receptors gabra5 and gabrg2, and the serotonin receptor htr1b were the three SDMN/stress-related genes with the highest number of direct connections to odorant and vomeronasal receptor genes (based on ranked ‘degree’ scores obtained from the Network Analyzer tool in Cytoscape; Assenov et al., 2007; Shannon et al., 2003) (Table S6).

Female conspecific mate preference in X. birchmanni requires experience with social cues. Although both conspecific- and heterospecific-exposed females developed a preference, consistent with previous studies (Verzijden and Rosenthal, 2011; Verzijden et al., 2012b), females socially isolated from adults failed to develop a preference. Experience with adults therefore appears necessary for the development of mating preferences (Fig. 1A).

Social experience with heterospecifics also shapes personality (Fig. 1B,C): X. malinche-exposed individuals resembled their shy tutors along boldness measures. In contrast to mate choice, both conspecific-exposed females and no-exposure fish showed boldness measures more like the conspecific X. birchmanni. Female X. birchmanni reared with heterospecifics therefore express mating biases and personality traits characteristic of the adults in their environment, which may facilitate rampant hybridization between these sister species that overlap in geographical distribution (Schumer et al., 2017; Culumber et al., 2011; Rosenthal et al., 2003). The fact that hybridization occurred only recently between these species (Schumer et al., 2018) could reflect the rarity of social environments needed to result in hybridization. For example, the most common form of migration according to species in this system is likely downstream migration of juvenile X. malinche, as they are found in the most upstream portions of each river (Culumber et al., 2011), and flooding and other disturbances predominantly disperses juvenile poeciliids in the direction of the river flow (Chapman and Kramer, 1991). As such, the more likely social scenario among species is for a cohort of relatively few X. malinche juveniles to be exposed to a predominantly X. birchmanni social environment. In this scenario, juvenile X. malinche should be expected to develop a relative disdain for the common X. birchmanni (Cui et al., 2017), while X. birchmanni develop a reinforced preference for the more familiar conspecifics, as seen in this and a previous study (Verzijden and Rosenthal, 2011), thereby reinforcing reproductive isolation. However, this expectation is based on the findings from studies focusing on exclusive exposure to one species versus the other. Future studies should determine how varying ratios of a mixed social environment according to species affects the development of these socially relevant behaviors.

Cultural transmission of personality traits may be adaptive (Schuett et al., 2013, 2010; Keiser et al., 2015), especially in variable environments or environments in which genetic parents are frequently separated from offspring or provide little to no care upon birth, a common occurrence among poeciliids (Chapman and Kramer, 1991). In these scenarios, it would be more informative for individuals to learn from the current environment's stock of adults instead of having personality-related traits be genetically predetermined from the parents. Learning from adults, regardless of genetic relationship, could provide a mechanism for the development of adaptive and phenotypically plastic personality traits.

Relative to socially exposed subjects, socially isolated females showed wide variance in expression profiles, with relative upregulation of both sensory processing and behavior-related genes. In addition to downregulating sensory genes, heterospecific-exposed females showed downregulated immune genes. This is in accord with previous research (Louder et al., 2018), consistent both with a role in learning and memory (Marin and Kipnis, 2013) and anticipation of sexual threat (Lawniczak and Begun, 2004; Bailey et al., 2011; Immonen and Ritchie, 2011) (Figs 3 and 4).

Heterospecific exposure was associated with the downregulation of an entire suite of odorant-related genes relative to both socially isolated and conspecific-exposed females (Fig. 4B,D), as well as of MHC-1-like genes (Fig. 4A,C, Table S4), which have important functional roles in both immune function and social communication (Milinski et al., 2005; Roberts and Gosling, 2003; Winternitz et al., 2017), and have previously been shown to have a negative correlation with synaptic plasticity in certain brain regions (Marin and Kipnis, 2013; Syken et al., 2006), as seen in our results. These results are in stark contrast to those seen in the sister species, X. malinche, where only one odorant receptor was found to be differentially expressed between exposure groups (Cui et al., 2017). Overall, we found few genes that were differentially expressed between social treatments in both studies, although 12 of these 16 genes do have functional roles in memory consolidation, learning, neurocognitive function or addiction. This list is also enriched for genes associated with the lateral septal nucleus (ventral telencephalic nucleus in teleosts), a brain region serving as a functional connection between the mesolimbic reward system and the social behavior network (O'Connell and Hofmann, 2011), though it is worth noting that these are based on the constitutive expression patterns of homologous human genes (Table S3). This low abundance of overlapping genes also suggests that the detection and learning of socially relevant signals in X. birchmanni and X. malinche may largely rely on different molecular pathways, which could explain the drastic differences in behavioral responses regarding social learning.

Whereas odorant-related genes were downregulated, the 451 genes upregulated in heterospecific-exposed females largely pertained to learning and memory consolidation processes (Bruel-Jungerman et al., 2007) such as neurogenesis and synaptic plasticity (Fig. 3B). Contrary to the overall pattern, a subset of SDMN and stress-related genes within the SDMN/stress module, most notably nlgn2b and npy8ar, showed strong co-expression links with genes within the vision and olfaction modules (Fig. 5). Inhibiting NLGN-2 protein increases boldness in mice (Liang et al., 2015; Kohl et al., 2015), consistent with reduced nlgn2b expression and increased boldness in X. birchmanni-exposed females. Neuroligin expression across the brain also increases with strength of preference during immediate mate choice in poeciliids (Wong and Cummings, 2014). Therefore, this family of genes may play an important role in the behavioral coupling of mate choice and personality-related traits. Furthermore, odorant receptor genes showed strong connections with stress-related genes with known anxiety-related functions [htr1b, gabrg2 and gabra5 (Lo Iacono et al., 2016; Backström and Winberg, 2017; Purdy et al., 1991), Fig. 6]. These results highlight a direct neurogenomic connection between sensory detection and the expression of social behaviors in response to differences in chronic social exposure.

The wide variance in NO-EXP transcriptomic profiles across multivariate space (Fig. 2A) may reflect an absence of structure that social exposure normally provides during neural development (Ten Cate, 1994; Marler, 1991), with differences in exposure altering the neurogenomic trajectory of female X. birchmanni brains. In zebra finches, social exposure during critical time periods results in a neuronal ‘shedding’ process in key brain regions, where only those neurons that are activated in response to a given social stimulus are retained, and inactive neurons are lost (Bischof, 2003; Changeux and Mikoshiba, 1978). Although differences in brain development between socially exposed versus isolated individuals have been described, few predictions can be made on how differences across exposure types can be expected to affect neural development (but see Louder et al., 2018). Our findings suggest that exposure to an adult model largely results in the relative prioritization of key molecular pathways and that the pathways that become prioritized are determined by the genetic relationship between the focal individual and the model. Hence, a lack of adult exposure altogether results in no change in the baseline response of these socially sensitive molecular pathways. For example, socially isolated females largely showed intermediate expression values of SDMN/stress-related genes compared with females of either exposure type (Fig. 4C). This intermediate expression across genes differentially expressed between B-EXP and M-EXP females would consequentially result in greater overlap between NO-EXP females and females from either exposure type than between the exposure types (B-EXP versus M-EXP) themselves. It would be interesting to determine the extent to which this ‘lack of response’ can be rescued via consistent exposure to an adult model later in life. Female mating preferences in this system can be reversed upon short-term exposure (Verzijden et al., 2012b), but does this behavioral change in preference require regulatory changes to the same pathways influenced by long-term differences in social exposure?

As with learned preferences for conspecifics in birdsong (Braaten and Reynolds, 1999; Doupe and Kuhl, 1999; Campbell and Hauber, 2009), X. birchmanni may be tuned to respond to sensory cues from conspecifics. Lack of exposure to key adult cues, or exposure to a novel phenotype, could result in a ‘rewiring’ process requiring a relative change in the regulation of key pathways related to social behavior (Louder et al., 2018). This effect of the social environment on baseline expression levels of social behavior genes could have important implications on the immediate neural mechanisms through which mating decisions are made (Bloch et al., 2018; Cummings, 2015), as several show baseline differences between social exposure types in this study (Table S4).

We found that female X. birchmanni develop both mating preferences and personality-related traits typical of the species to which they were exposed. Shared cultural transmission of personality and homotypic mating preferences could function as a ‘magic trait’ (Servedio et al., 2011), coupling assortative mating with ecological divergence. Our results suggest that this coupling is due not only to correlations between sexual signals and personality traits in adult fish that serve as models, but also to neural mechanisms that link an individual's perceived social environment to their preferences and personality. In this study, we described long-term changes in transcriptomic profiles of the female brain and sensory periphery in response to differences in the social environment, though we acknowledge that our results fail to identify site-specific changes in gene expression within brain regions. Although conspecific exposure resulted in an upregulation of immune response and sensory detection genes, heterospecific exposure led to an upregulation of genes related to neurogenesis and synaptic plasticity, suggesting an intriguing split in the neurogenomic trajectory of individuals with differing social pasts. Mirroring the behavioral link we observed between mating preference and shyness–boldness, we found that certain stress-coping genes and SDMN pathway genes were found within the same module and showed strong co-expression with genes related to sensory detection. Specifically, neuroligins may shape both female mating preferences and personality, and we found nlgn2b to play a central role in distinguishing between the expression profiles of conspecific- versus heterospecific-exposed females. Neuroligins are thus promising candidates for future research aiming to determine how ‘magic traits’ (Servedio et al., 2011) are developed and shaped by the social environment. Furthermore, our results show that social exposure has long-term consequences on the baseline expression profiles of molecular pathways previously shown to have functional roles in the immediate mate choice process (Bloch et al., 2018; Cummings, 2015; Wong et al., 2012; Wong and Cummings, 2014), thus highlighting an important interplay between long- and short-term social experience in shaping socially relevant behaviors. Specifically, immediate gene expression and behavioral responses to social stimuli may strongly depend on the baseline transcriptomic profile of an individual, which is a product of their social past.

We thank Angela Freeman, Molly Schumer, Rongfeng Cui, Daniel Powell, Gastón Jofre and Christopher Holland for comments on the experimental design and manuscript. We thank Molly Cummings for discussion of topics directly related to the manuscript. We also thank Geneva Kernaghan and Abby Glueck for aid in conducting behavior trials and video scoring.