<figure><img src="image_1.png" /></figure>

HAL
open science

# Kin recognition: Neurogenomic response to mate choice and sib mating avoidance in a parasitic wasp

Aurore Gallot, Sandrine Sauzet, Emmanuel Desouhant

► To cite this version:

Aurore Gallot, Sandrine Sauzet, Emmanuel Desouhant. Kin recognition: Neurogenomic response to mate choice and sib mating avoidance in a parasitic wasp. PLoS ONE, 2020, 15 (10), pp.e0241128. 10.1371/journal.pone.0241128. hal-02984453

HAL Id: hal-02984453

https://cnrs.hal.science/hal-02984453

Submitted on 1 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
RESEARCH ARTICLE

# Kin recognition: Neurogenomic response to mate choice and sib mating avoidance in a parasitic wasp

Aurore Gallot<sup>1*</sup>, Sandrine Sauzet<sup>1,2</sup>, Emmanuel Desouhant<sup>1</sup>

1 Laboratoire de Biométrie et Biologie Evolutive, CNRS, Université Lyon 1, Université de Lyon, UMR 5558, Villeurbanne, France, 2 Institut de Génétique Humaine, CNRS-Université de Montpellier, UMR 9002, Biology of Repetitive Sequences, Montpellier, France

* aurore.gallot@univ-lyon1.fr

<figure><img src="image_2.png" /></figure>

## Abstract

Sib mating increases homozygosity, which therefore increases the risk of inbreeding depression. Selective pressures have favoured the evolution of kin recognition and avoidance of sib mating in numerous species, including the parasitoid wasp *Venturia canescens*. We studied the female neurogenomic response associated with sib mating avoidance after females were exposed to courtship displays by i) unrelated males or ii) related males or iii) no courtship (controls). First, by comparing the transcriptional responses of females exposed to courtship displays to those exposed to controls, we saw a rapid and extensive transcriptional shift consistent with social environment. Second, by comparing the transcriptional responses of females exposed to courtship by related to those exposed to unrelated males, we characterized distinct and repeatable transcriptomic patterns that correlated with the relatedness of the courting male. Network analysis revealed 3 modules of specific 'sib-responsive' genes that were distinct from other 'courtship-responsive' modules. Therefore, specific neurogenomic states with characteristic brain transcriptomes associated with different behavioural responses affect sib mating avoidance behaviour.

<figure><img src="image_3.png" /></figure>

### OPEN ACCESS

**Citation:** Gallot A, Sauzet S, Desouhant E (2020) Kin recognition: Neurogenomic response to mate choice and sib mating avoidance in a parasitic wasp. PLoS ONE 15(10): e0241128. https://doi.org/10.1371/journal.pone.0241128

**Editor:** William J. Etges, University of Arkansas, UNITED STATES

**Received:** January 17, 2020

**Accepted:** October 8, 2020

**Published:** October 26, 2020

**Peer Review History:** PLOS recognizes the benefits of transparency in the peer review process; therefore, we enable the publication of all of the content of peer review and author responses alongside final, published articles. The editorial history of this article is available here: https://doi.org/10.1371/journal.pone.0241128

**Copyright:** © 2020 Gallot et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

**Data Availability Statement:** Raw sequencing data, head transcriptome data, the expression raw count table per gene (HTSEQ output) and associated metadata are available in the NCBI Gene Expression Omnibus (GEO) repository and

## Introduction

Mate choice often depends on individual attractiveness. In many species, females discriminate between competing males [1]. Female mate choice sometimes occurs after elaborate courtship displays, and males are selected according to their highly heritable qualities ('good genes' hypothesis) and/or their genetic compatibility with the female ('compatible genes' hypothesis) [2]. The 'good genes' hypothesis predicts that females favour reproduction with males carrying traits that are honest indicators of good genes or as a result of sensory bias [3, 4], hence obtaining genetic benefits for their offspring [5]. The 'compatible genes' hypothesis predicts that offspring fitness is correlated with genetic dissimilarity between mating partners, where each female prefers males that possess genes compatible with their own genotypes. Sib mating avoidance fits into this second hypothesis [6, 7]; in this case, kin recognition is crucial for enabling females to discriminate between related and unrelated mates. Sib mating avoidance
accessible through accession number GSE124129. Normalized expression dataset, matrix reporting the experimental design, and the script used to analyzed dataset are available as Supporting Information. A full list of differentially expressed genes and associated statistics is available in S2 and S3 Tables. Additional data may be requested from the authors.

**Funding:** This work was funded by the Agence Nationale de la Recherche (https://anr.fr/en/) JCJC AVOIDINBRED (ANR-17-CE02-0004-01) attributed to AG. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

**Competing interests:** The authors have declared that no competing interests exist.

has been extensively documented in numerous animal species (including species that are social and solitary) and reduces the risk of inbreeding depression [8].

In colonial marine invertebrates, contact between 2 conspecific colonies may trigger fusion if the 2 colonies are identical or related, or aggressive rejection if the 2 colonies are not related. Molecular bases for kin recognition involve several highly polymorphic loci that have been well characterized in 2 species: the cnidaria *Hydractinia symbiolongicarpus* (reviewed in [9]) and the chordata *Botryllus schosseri* (reviewed in [10]). In terrestrial animals, most of the presently identified mechanisms usually rely on individual-specific olfactory cues and sensory organs able to discriminate these cues from one another. In vertebrates, major histocompatibility complex genes confer individual olfactory identities [11–14]. These olfactory cues bind to receptors often located on the vomeronasal organ neurons and are used by females to evaluate mate relatedness [15–19]. In insects, numerous species discriminate kin by using odour differences inherent to divergence in the cuticular hydrocarbon composition [20–23]. Odour molecules bind to specific sensory receptors located on sensilla basiconica, where dozens of olfactory sensory neurons project to the antennal lobes [24–26], where olfactory second-order neurons project to the lateral horn and mushroom bodies [27, 28]. To the best of our knowledge, no genome-wide study has addressed how kin recognition signals are processed in the context of mate choice and sib mating avoidance. Transcriptomic studies provided recent insight into female mating decisions [29]. Coordinated changes in the expression of many genes in female brains, *i.e.*, a neurogenomic response, have been identified following courtship in Poeciliidae fishes [30–32]. This response depends on male attractiveness and is in accordance with female preferences.

Hymenopterans (including bees, ants, wasps and sawflies) are relevant models to explore the biology of mate choice, especially regarding sib mating avoidance and its underlying biochemical mechanisms. These insects share an ancestral and unusual sex determination system called single-locus complementary sex determination (sl-CSD), where sex is determined by heterozygosity at the sl-SCD locus. Heterozygous individuals develop as females, hemizygous individuals develop as males, but homozygous individuals develop into diploid unviable or sterile males [33–36]. When sib mating occurs, the risk of genetic incompatibility between partners sharing the same allele at the sl-CSD locus is one in two; in that case, half of the progeny will be homozygous diploid sterile males (with notable exceptions [37, 38]). Therefore, the sl-CSD has an additional immediate effect on inbreeding. Consequently, Hymenoptera are much more exposed to inbreeding depression than any other diploid species [39]. In this context, congruent with the compatible gene hypothesis, one expects that selection will favour females that prefer males with a distinct allele at the sl-CSD locus, thus excluding sibs. In the parasitoid wasp *Venturia canescens*, which has sl-CSD [40], females only mate once [41], making mate choice particularly decisive. Indeed, in this species, females are able to discriminate kin and non-kin during male courtship based on olfactory-mediated cues [42]. In laboratory conditions, the proportion of successful mate when a single female was in presence of 2 unrelated males was higher (79% of success) than in presence of related males (*i.e.* with 2 brothers) (41% of success) [42]. Furthermore, the latency to mate increased when the courting male was related [42]. Experiments involving choice between 2 males, one brother and one unrelated present in a same area have shown that female mates indifferently with the 2 males [42]. Hence, the two distinct social contexts (presence of a conspecific with and without relation) provoke two contrasting behavioural responses that correspond to the definition of behavioural states, *i.e.*, the performance of a distinct and quantifiable behaviour for a measurable period of time [43].

We tested whether these contrasting behavioural responses to relatedness were correlated to transcriptomic changes in *V. canescens* females. We designed a behavioural experiment
where females were exposed for 10 minutes to an unrelated male or a related male (a brother) or were without any social contact. Then, we used RNA sequencing to compare gene expression profiles from the entire heads of these females. This timing coincided with a period of active courtship, where female wasps perceived and evaluated males, experienced change in receptivity and decided whether or not to mate [44]. Previous studies demonstrated that this temporal frame is adequate to detect early transcriptional changes following mate exposure in different species [30, 45, 46]. We focused on the whole head, where sensory signals are processed and integrated to mediate complex decision making, including mating decisions.

We hypothesized that information processing during male courtship display and mate choice were related to gene expression changes. First, we examined if courtship perception was mediated by changes in gene expression. In this case, we expected the transcriptomic profiles of all courted females would be analogous, whatever the relatedness of the courting male. Then, we assessed whether the differences in female receptivity were associated with gene expression changes. Under this hypothesis, the expectation was that non-receptive females, *i. e.*, solitary females and females courted by related males, would have similar transcriptomic profiles.

# Materials and methods

## Biological model

*Venturia canescens* G. (Hymenoptera: Ichneumonidae) is a solitary parasitic wasp; a maximum of one adult emerges from one host, regardless of the number of parasitoid eggs initially laid in the host. In the wild, *V. canescens* females parasitize pyralid moth larvae feeding on dried fruits, such as figs, carobs, almonds, dates or loquats [47]. Our knowledge of how mating partners encounter each other in the field is largely incomplete as a consequence of the small size of the species, which renders observations difficult. Virgin *V. canescens* females emit chemicals that, in combination with host kairomones, attract males [42]. In turn, males do not attract virgin females at a distance [41]. *V. canescens* females are monandrous and only mate once in their lifetime [41]. Conversely, males can mate more than once. These observations led us to consider that females choose their mate, and we thus focused our study on female transcriptomic responses to courtship. Moreover, females recognize sibs on the basis of a chemical signature emitted by males and can avoid mating with their brothers through kin recognition [42]. This kin recognition ability has also been robustly established in this species in the context of host choice, since females prefer to lay their eggs in hosts already parasitized by conspecifics rather than in hosts already parasitized by themselves [48]. In addition, males also have the ability to distinguish between non-sib and sib females by using chemical markers emitted by females [49].

## Behavioural manipulations

Wild parasitic wasps were collected via a large field sampling during the summer of 2014 in orchards in southern France (Valence, N 44°58'21" E 4°55'39"). Wasps were maintained on the host *Ephestia kuehniella* (Lepidoptera: Pyralidae) in a climatic chamber (25°C, 60% humidity, DL 12:12) and fed with a 50/50 honey/water solution. Families were produced as described by Metzger et al. (2010). A total of 90 newly emerged virgin females were individually placed in a box and randomly submitted to one of the three following conditions for a maximum duration of 10 minutes: 1) a compatible male (*i.e.*, unrelated male), 2) an incompatible male (*i. e.*, a brother), or 3) isolated (controls). The experimental design has been conceived in order to minimize the influence of circadian rhythm and genetic background on the results. Newly emerged females were isolated every morning and numbered. A random draw was then made
to establish the order of passage of the females and the treatment assigned. Behavioural experi-
ments took place in the afternoon between 1:00 and 4:00 p.m. for all females that were cap-
tured in the morning in the order established by the random draw. A maximum of 2 females
from each family have been kept daily, and randomized in one condition. For each given con-
dition, the ten females belonged to different families to avoid an effect of genetic homogeneity
on transcriptomic results. New families have been produced for every biological replicate. The
10 minutes period coincides with an active male courtship behaviour [44]. For each female
under conditions 1 and 2, one sequence of active courtship parade was observed within the 10
minutes following the introduction of the male to the box for all pairs. Whenever this behav-
iour sequence was not observed, the pair was eliminated from the study. After this time, or
within 3 seconds of mounting, the female and male were separated. For each condition,
females were caught in a small tube (3 x 7 cm) and were immediately anaesthetized with CO2
and decapitated with a scalpel. Antennae were removed from the head to focus on the cerebral
transcriptome of the decision-making centre rather than the sensory centre. The head is a het-
erogeneous tissue divided into many functionally distinct regions that could affect inferences
about the functional significance of gene expression patterns [50]. However, the study of more
specific brain regions remains challenging given the small size of the wasps. Collected heads
were immediately flash-frozen and stored at -80°C until RNA extraction. A total of 10 individ-
ual heads were collected per experimental condition and pooled to collect higher RNA
amounts and to average the expression state across individuals, thus mitigating individual vari-
ability during transcriptome comparison [51]. The experiment was sequentially repeated 3
times to obtain biological replicates.

RNA extraction, library preparation and sequencing

Total RNA extraction was performed in one batch on the 9 biological samples collected (3 experimental conditions x 3 biological replicates). Mechanical cell lysis was performed with metallic beads added to frozen microtubes and shaken with a TissueLyser (Qiagen, 25 hertz 45 seconds). RNA was extracted using Qiazol and an RNeasy mini kit according to the manufacturer's protocol (Qiagen), with optional DNase treatment (Thermo Fisher Scientific). The quality and quantity of total RNAs were assessed on a Bioanalyzer (Agilent) and a Qubit fluorometer (Life Technologies). Polyadenylated RNAs were enriched from 1 µg of high-quality total RNAs with oligo-dT magnetic beads, fragmented and converted to cDNAs (Illumina Tru-Seq Stranded mRNA Library Prep kit). After adding an A to the 3' end of each cDNA, adapters were ligated, and fragments were amplified by PCR to generate DNA colonies. Each library was labelled, multiplexed and pooled for sequencing on a HiSeq 2500 Illumina sequencer (Fasteris, Switzerland), with a paired-end protocol (2x125 bp). For each of the 9 biological samples, a second library was prepared independently, sequenced and considered a technical replicate. Overall, a total of 18 libraries were sequenced (3 experimental conditions x 3 biological replicates x 2 technical replicates). A minimum of 18 million paired-reads were obtained per library, representing a final dataset of more than 535 million paired-reads.

RNA-seq quality control, mapping, transcriptome assembly and
annotation

The *V. canescens* genome has been sequenced and is currently available (http://bipaa.genouest.org/sp/venturia_canescens/, v.1.0) but thus far, RNA-seq has not been performed on heads [52, 53]. To include all loci detected in the head transcriptome that might be missing from the current annotation, we constructed a *de novo* transcriptome assembly using the genome reference with the TopHat-Cufflinks pipeline [54] (v2.2.1). The global quality of sequences was
assessed with FastQC-0.10.1. Low-quality bases at the ends of reads were trimmed using
Trimmomatic-0.36. Adaptor-containing reads, low quality reads (scores phred <30) and N-
containing reads were filtered out. Curated reads were then processed through the TopHat-
Cufflinks pipeline. First, paired-reads of each sample were aligned to the *V. canescens*
genome. TopHat used the Bowtie algorithm to align reads to the genome. Unmapped reads
were cut into segments that aligned far apart from one another (between 100 bp to several
hundred kb) to predict potential intron/exon structures. Next, an index of hypothetical splice
sites was built without any prior information. Each site had to be confirmed by several read
segments consistently showing the same alignment pattern. TopHat was used with each of
the 18 libraries as input and with the following parameters: stranded libraries (first strand), 5
mismatches/indels allowed, and report the best alignment possible for each read. After map-
ping, resulting alignment files were provided to Cufflinks, which generated a transcriptome
assembly for each sample annotated into genes, transcripts and isoforms. Then, the 18 tran-
scriptomes were merged into one master transcriptome with Cuffmerge. This final transcrip-
tome assembly contained a total of 16,752 genes and provided a uniform basis for calculating
gene expression in each condition. The transcriptome was then annotated to gain insight into
the functions of the transcripts and the proteins they encode. Sequence similarity was
searched for each of the 16,752 predicted genes by comparing the six-frame translation puta-
tive products of the nucleotide sequences using BLASTX (v2.2.29+) against the NCBI non-
redundant protein database. Transcriptome analysis was completed using gene ontology
(GO) annotation, which associated genes with functions in 3 categories: molecular function
(molecular activities of gene products), cellular component (where gene products are active)
and biological process (pathways and larger processes made up of the activities of multiple
gene products). This classification enabled functional interpretation of a large group of genes
via enrichment analysis.

Sample-based clustering and PCA

To obtain an initial overview of gene expression patterns across samples, multivariate analyses
were performed on the 18 transcriptomic libraries, including 3 biological replicates for each of
the 3 conditions representing a total of 9 samples, each one represented by 2 technical repli-
cates. Those 18 transcriptomes were represented by isolated females (6), females courted by
unrelated males (6) or females courted by brothers (6). A raw count table was obtained by
using HTSeq [55] (v0.5.4p1). The gene model produced with Cuffmerge was combined with
the 18 mapping files previously obtained with TopHat for each gene in each library. This final
dataset was exported to R (v3.4.3) [56] for downstream statistical analysis with DESeq2 [57]
(v1.30.0). Counts were normalized with the variance stabilizing transformation method (VST),
which produced transformed data on the log2 scale and normalized data with respect to library
size. The overall variation of expression levels among samples was evaluated with a two-
dimensional PCA and with hierarchical clustering, both based on the expression of the 500
genes with the highest variance across libraries. A hierarchical clustering dendrogram was pro-
duced based on the sample-to-sample Euclidian distance matrix to obtain an overview of the
similarities and dissimilarities between samples. Uncertainty in hierarchical clustering was
assessed with Pvalclust [58] (v2.0-0) using multiscale bootstrap resampling with an approxi-
mately unbiased P-value to measure statistical support for each cluster (1,000,000 bootstraps;
average agglomerative method; correlation method distance). Two technical replicates from
one biological sample of females courted by brothers were obvious outliers in the PCA and
hierarchical clustering, most likely due to a problem that occurred during sample conditioning
and were thus excluded from the analysis (S1 Fig).
Differential expression, co-expression network and functional enrichment
analysis

All pairs of technical replicates were merged before proceeding to differential expression anal-
ysis and network analysis as recommended by Love et al. (2014), keeping a total of 8 samples
(3 isolated females, 3 females courted by unrelated males and 2 females courted by brothers).
To identify genes with different expression patterns across conditions, we performed pairwise
comparisons between 1) females courted by related and unrelated males to controls (i.e., iso-
lated females); and 2) females courted by unrelated males to females courted by related males.
Differential expression was tested by using negative binomial generalized linear models imple-
mented in the program DESeq2. We tested for differential expression of all transcripts with an
average level of expression superior to 10 reads per library (n = 14,034). After normalization,
within-group variability (i.e., the variability between biological replicates) was modelled for
each gene by the dispersion parameter, which describes the variance of counts by sharing
information across genes, assuming that genes of similar average expression strength have a
similar dispersion. Among the 14,034 genes tested, only 7 presented an outlier status (i.e., an
inconsistent expression pattern across the dataset) and were thus excluded from the test, keep-
ing 14,027 transcript in the final reference transcriptome. For the subset of genes that passed
the filtering test, a Wald test P-value was calculated and finally adjusted for multiple testing
[59]. A gene was considered differentially expressed when the false discovery rate (FDR)-
adjusted P-value was less than 0.01. We did not apply any log fold change threshold.

We performed a weighted gene co-expression network analysis (WGCNA [60], v1.63) to identify subgroups of genes that shared common expression patterns across the experimental conditions and potentially drove the differences in mate choice. WGCNA is a data reduction technique that regroups genes with similar expression patterns into modules of co-expressed genes and tests the correlations between modules and traits. First, log2-transformed and VST-normalized counts were used to construct a gene co-expression network with the blockwise-Modules function. A correlation matrix was computed for the genes, and the correlations were weighted using a power function $\beta$. Then, genes sharing similar patterns of variation across conditions were regrouped using hierarchical clustering and a dynamic tree-cutting algorithm to define modules of co-expressed genes. For our analysis, the parameters used were as follows: maximum block size = 15,000 genes; power ($\beta$) = 10, minimum module size = 30 genes. The remaining parameters were kept at the default settings. A colour name was assigned to each module, and biologically interesting modules were identified by correlating a summary profile for each module to external experimental conditions; P-values <0.05 were considered significant and numbered from 1 to 11. Finally, the potential 'hub genes' in every significant module were identified. So-called hub genes may influence the expression of other genes in their module and may be causal factors for a trait of interest. Such hub genes are potentially biologically relevant by driving phenotypic variations [61, 62]. The identification of hub genes relies on both connectivity with other genes from the module and the correlation to the trait. Accordingly, genes were ranked according to their module membership values in each module. The top 5 hub genes of every module were annotated, and their expression pattern, which was representative of the module they belonged to, was detailed.

Functional characterization of gene sets (i.e., genes with differential expression, or modules
of co-expressed genes) was analysed using enrichment analysis and gene ontology annotations.
An enrichment test was performed on test sets compared to the full transcriptome with Blas-
t2GO [63] (v5.0). The proportion of genes associated with GO terms was compared between
the test set and the transcriptome (14,027 transcripts) with a unilateral Fisher's exact test (one-
sided), P-values < 0.01 were considered significantly enriched.
# Results

## The transcriptomic response to courtship is consistent with the mate preference of females

We first determined whether the male courtship display provoked a transcriptomic response in females. To do so, we built the first *V. canescens* head transcriptome by sequencing more than 500 million reads using a genome guide assembly (18 to 55 million paired-reads from 18 libraries, S1 Table) and used it as a reference to compare the different head transcriptomes of females. After quality filtering a mean of 98.4% of paired-reads were kept, on which 70.8% were successfully mapped to the genome (representing a mean of 20 millions per sample, for a total of 363 million paired-reads, S1 Table). Such values corresponded to the high-quality standards observed in other Hymenopteran species with an annotated genome [64]. The transcriptome constructed with these sequences encompassed 16,752 genes. Overall, 76% of predicted genes get a blast hit (12,740) while 4,012 sequences get no hit. Among genes matching with blast, 89.4% presented their best hit with an insect sequence, of which 84.4% match more specifically to a hymenopteran insect sequences. Finally, most of the genes were uncharacterized, since only 33.4% of the predicted genes (5,589) were associated with at least one Gene Ontology (GO) functional annotation.

Based on this unique reference, head transcriptomes from females 1) courted by an unrelated male, 2) courted by a related male, and 3) isolated (controls) were first analysed without a priori knowledge. Principal component analysis (PCA, Fig 1A) revealed that transcriptomes shared higher similarity within a given experimental condition than between different experimental conditions. PCA defined three consistent clusters according to the social environment proposed: the control group where females were kept in isolation, the group of females courted

<figure><img src="image_4.png" /><figcaption>**Fig 1.** Multivariate analysis of the 16 RNA-seq libraries based on the gene expression profiles of the 500 genes with the highest variance across samples. The samples are clustered by the experimental conditions that females were exposed to. Experimental conditions are represented by the three different colours (grey: control; orange: females courted by related males; blue: courted by unrelated males), shapes indicate the three biological replicates, and the filling of the shape shows the 2 technical replicates. (A) Two-dimensional principal component analysis (PCA) plot. (B) Hierarchical clustering dendrogram based on sample-to-sample distances. Statistical support is indicated by an approximately unbiased *P*-value with one million multiscale bootstrap resamplings (all bootstrap values were > 50%, those < 80% were not shown for clarity).</figcaption></figure>

https://doi.org/10.1371/journal.pone.0241128.g001
by unrelated males, and the group of females courted by related males. Principal component 1 separated females courted by related males from the two other conditions and explained the largest fraction of variance in gene expression (46%). Principal component 2 separated females courted by unrelated males from the two other conditions and accounted for 17% of the variance in gene expression (Fig 1A). These results were supported by hierarchical clustering analysis of sample-to-sample distance, which established that samples were clustered by the experimental conditions that they were exposed to (Fig 1B). Bootstrap resampling provided strong statistical support for this result: the control group (isolated females) formed a cluster with 100% support (Fig 1B), while the group of females courted by unrelated males constituted a cluster with 100% support. Finally, females courted by related males constituted a group with 95% support (Fig 1B). Together, these results showed that i) male courtship had an effect on the female transcriptome, and ii) being courted by a related or an unrelated male provoked two distinct transcriptomic responses.

## Functional characterization of differentially expressed genes following courtship

The transcriptomes of females courted by related and unrelated males were compared to those of isolated females to identify the female neurogenomic response to courtship. We identified 1,001 differentially expressed genes (DEGs), representing 7.1% of the 14,027 gene sets tested (Fig 2A, S2 Table). Among the 1,001 DEGs, 463 had higher expression in isolated females (3.3% of total transcriptome), gene ontology enrichment analysis reveals that this set of gene was enriched in DNA-binding Transcription Factor Activity (full list in S2 Table). In contrast, 538 DEGs were overexpressed in courted females (3.8% of total transcriptome), gene ontology enrichment analysis reveals that this set of gene was enriched in Reproductive Behaviour (full list in S2 Table).

Among the DEGs, we noticed three neuropeptides associated with female receptivity, namely, *Capability*, *PDF* and *SIFamide*, were all downregulated in courted females (Fig 2). The orthologue to *muscleblind*, required for normal regulation of female sexual receptivity [65], was also detected among the genes overexpressed in courted females. Notably, we also reported two genes associated with dopamine, *ebony* (downregulated in courted females), *DAT* (*dopamine transporter*, overexpressed in courted females), and one gene associated with serotonin (*5-HT1A*, *serotonin receptor 1A*, overexpressed in courted females) (Fig 2).

## Relatedness of the courting male influences the female head transcriptome

We compared females courted by unrelated males to those courted by related males to determine whether the relatedness to the courting male had an impact on the female transcriptomic response. By comparing these two groups, we identified 831 DEGs representing 5.9% of the tested genes (Fig 3A, S3 Table). Among these DEGs, 481 genes presented higher expression in females courted by related males (3.4%), gene ontology enrichment analysis reveals that this set of gene was enriched including ATP metabolism and Ribosome (full list in S3 Table). Moreover, 350 DEGs were overexpressed in females courted by unrelated males (2.5%). Enrichment analysis revealed 22 associated GO terms including Reproductive Behaviour and Male Mating Behaviour (full list in S3 Table). These 2 GO terms both refer to the same 2 genes, which belong to the *yellow-major royal jelly protein* family (Fig 3B, S2 Fig), an insect gene family notably associated with behaviour and caste specification [66]. Interestingly, we also noticed the regulation of the transcription factor *PAX6*, which is required for normal brain structure and function, notably locomotory behaviour [67] (Fig 3B).
A

<figure><img src="image_5.png" /></figure>

B

<figure><img src="image_6.png" /></figure>
**Fig 2.** Transcriptomic response to courtship: 7.1% of the total transcriptome (grey) was regulated in response to courtship (1,001 DEGs). (A) Comparison of females courted by related and unrelated males to controls revealed 538 genes were overexpressed in courted females (negative fold change values) and 463 genes were overexpressed in isolated females (positive fold change values). (B) Boxplot showing significant expression changes following courtship in highlighted genes mentioned in the text. The X axes indicate gene names, the Y axes show normalized counts after log transformation, and the boxplot whiskers show the range of read counts between biological replicates. PDF, pigment-dispersing factor; DAT, Dopamine transporter; 5-HT1A, Serotonin transporter 1A.

https://doi.org/10.1371/journal.pone.0241128.g002

## Co-expression network and characterization of ‘courtship-responsive’ modules and ‘sib-responsive’ modules

To better characterize the variations in gene expression according to the experimental conditions experienced by females, we applied a co-expression network analysis using WGCNA on the 14,027 genes that passed through the expression filter. The gene co-expression network grouped genes that shared a similar expression pattern across different experimental conditions into modules.

Overall, the 14,027 genes were organized into 50 modules of highly correlated genes symbolized by a colour, with sizes varying from 32 to 2,675 genes (Fig 4, Table 1). Among the 50 modules defined by the cluster analysis, 11 modules (numbered from 1 to 11) had significant correlations with at least one of the experimental conditions experienced by the females (Table 1, Fig 4B, S3 Fig); these were considered biologically relevant and were further analysed.

First, three modules were associated with courtship display by an unrelated male (1, 2 and 3), which grouped 1,239 genes. Then, three modules were associated with sib-responsive genes (4, 5 and 6; S3 Fig; S4 Table), which grouped 2,824 genes responding only in the presence of related males. Next, two modules were associated with response to courtship, whatever the degree of relatedness of the courting males (modules 7 and 8; S3 Fig; S4 Table), which grouped 2,780 genes. Finally, three modules were associated with a gene expression pattern peculiar to each of the 3 social environments (9, 10 and 11; S3 Fig; S4 Table).

## Discussion

In the current study, we characterized the female neurogenomic response associated with sib mating avoidance by identifying remarkable differences in the head transcriptome triggered by courtship display that differed according to the relatedness of the courting male. In *V. canescens* females, mate relatedness influences female sexual receptivity and is estimated during male courtship displays through chemical cues [42]. Unrelated males induce sexual receptivity in females, whereas related males induce weak sexual receptivity. Hence, sib mating avoidance can be considered a behavioural state, similar to other transitory behaviours such as aggressiveness [68] or singing [69]. We showed that sib mating avoidance is associated with distinct and repeatable cerebral transcriptomic patterns involving a significant part of the transcriptome (>5%). Such results fit the definition of neurogenomic state, i.e., a distinct and repeatable pattern of gene expression in the brain revealed by contrasting brain transcriptomes of individuals across different behavioural states. Despite the quite low number of biological replicates, the highly contrasted transcriptomes observed in the different social contexts suggest that sib mating avoidance behaviour could be considered a neurogenomic state. This research paves the way for further study on neurogenomic effects of sib mating avoidance in many species where such behaviours have been described and, thus, may contribute to the understanding of the molecular mechanisms underlying the evolution of avoiding consanguinity.

We measured major transcriptomic modifications occurring in the female head following a courtship display. Females exhibited transcriptomic pattern changes following an encounter
<figure><img src="image_7.png" /><figcaption>**Fig 3.** Influence of relatedness of courting males on females' transcriptome expression: 5.9% of the total transcriptome (grey) was regulated according to the relatedness of the courting male. (A) Comparison of females courted by related males to females courted by unrelated males showed 481 genes overexpressed in females courted by related males (positive fold change values) and 350 genes overexpressed in females courted by unrelated males (negative fold change values). (B) Boxplot showing significant expression changes according to the relatedness of the courting male for highlighted genes mentioned in the text. The X axes indicate gene names, Y axes show normalized counts after log transformation, and boxplot whiskers show the range of read counts between biological replicates. *Yellow—mrjp*: yellow major royal jelly proteins; *PAX6*: paired box 6.</figcaption></figure>

https://doi.org/10.1371/journal.pone.0241128.g003
A

<figure><img src="image_8.png" /></figure>

B

<figure><img src="image_9.png" /></figure>
**Fig 4.** Gene co-expression network of the full transcriptome (14,027 genes) in response to courtship and according to relatedness with the courting male. (A) Clustering dendrogram of 14,027 genes (top side), with assignation to 50 modules of co-expressed genes represented by colors (bottom side). B) Summary plot of modules dendrogram and relationship with social environment experienced. The right panel shows the dendrogram of modules from gene network. The heatmap (left panel) illustrates pairwise correlation between modules of co-expressed genes and social environment: red denotes high positive correlation, white indicates no correlation, while blue indicates high negative correlation. The 11 modules significantly correlated with at least one of the social environment experimented by females were numbered from 1 to 11. Social environments were labelled with symbols: §; control (isolated females), $; females courted by related males, and &; females courted by unrelated males.

https://doi.org/10.1371/journal.pone.0241128.g004

with a partner, regardless of the relatedness of the courting male. This suggests that transcriptomic shifts immediately arise following the male display, even before the eventual copulation. These transcriptomic changes can be triggered by the presence of a male, mate evaluation, or from a social encounter. With the current experimental protocol, it is not possible to distinguish the origin of these changes. The addition of two other controls such as females in the presence of i) an unrelated female or ii) a related female would allow further specification of the impact of social environment on neurogenomic responses. Overall, 7.1% of the total transcriptome was differentially expressed at most within ten minutes after the courtship started (1,001 DEGs). Such a neurogenomic response on the time scale of minutes following environmental change is mediated by immediate early genes [70]. For such genes, near-instantaneous transcription is allowed by the presence of RNA polymerase II, which stalls in the promoter regions of these genes [71].

A recent transcriptomic study conducted on female mate preference in the guppy (*Poecilia reticulata*) showed that the presence of a potential mate induced changes in the brain transcriptome after only 10 minutes of exposure [30]. In insects, few studies have attempted to identify transcriptional changes associated with courtship displays; to the best of our knowledge, all have focused on *Drosophila* [45, 46]. Immonen and Ritchie submitted *D. melanogaster* females to a courtship song diffused by a speaker for 15 minutes and found only 41 DEGs

<table><caption>Table 1. Significant correlations between 11 modules and courtship experienced by females.</caption><thead><tr><th>Module name</th><th>Module size</th><th>representative enriched GO terms</th><th>Controls (isolated females)</th><th>Females courted by related males</th><th>Females courted by unrelated males</th></tr></thead><tbody><tr><td>1 (purple)</td><td>227</td><td>molybdenum ion binding</td><td>NS</td><td>NS</td><td>-0.79 (0.02)</td></tr><tr><td>2 (green)</td><td>935</td><td>nucleotide metabolic process</td><td>NS</td><td>NS</td><td>-0.73 (0.04)</td></tr><tr><td>3 (pale turquoise)</td><td>77</td><td>mitochondrion, respiration</td><td>NS</td><td>NS</td><td>-0.8 (0.02)</td></tr><tr><td>4 (blue)</td><td>2108</td><td>ribosome, sensory perception</td><td>NS</td><td>0.94 (5e-04)</td><td>NS</td></tr><tr><td>5 (red)</td><td>643</td><td>Methylation</td><td>NS</td><td>-0.76 (0.03)</td><td>NS</td></tr><tr><td>6 (dark magenta)</td><td>73</td><td>GDP-mannose metabolic process</td><td>NS</td><td>-0.74 (0.04)</td><td>NS</td></tr><tr><td>7 (turquoise)</td><td>2675</td><td>mating behaviour</td><td>0.91 (0.002)</td><td>NS</td><td>NS</td></tr><tr><td>8 (dark green)</td><td>105</td><td>cell surface receptor signaling pathway</td><td>-0.77 (0.03)</td><td>NS</td><td>NS</td></tr><tr><td>9 (floral white)</td><td>38</td><td>oxidoreductase activity</td><td>NS</td><td>-0.8 (0.02)</td><td>0.82 (0.01)</td></tr><tr><td>10 (brown)</td><td>1989</td><td>response to stimulus, protein kinase activity</td><td>-0.73 (0.04)</td><td>NS</td><td>0.81 (0.01)</td></tr><tr><td>11 (yellow)</td><td>1022</td><td>carbohydrate catabolic process, metal ion binding</td><td>-0.8 (0.02)</td><td>0.79 (0.02)</td><td>NS</td></tr></tbody></table>

NS: not significant

Correlation coefficients and *P*-values in parentheses are indicated when significant (*P*-value <0.05). The results of enrichment analyses are indicated in the GO terms column (full lists in S4 Table).

https://doi.org/10.1371/journal.pone.0241128.t001
(0.3% of the transcriptome) between courted and control females (in the presence of a male that was unable to perform a courtship display). Veltsos et al. (2017) submitted *D. pseudoobscura* females to male courtship and examined head transcriptomes immediately after mounting. They identified only 16 DEGs between courted and control females (virgin isolated females) (0.1% of the transcriptome). Together, these studies suggest that the female neurogenomic response to courtship only affects a very small set (>1%) of genes in Drosophila. Interestingly, in contrast with those studies, we found that the neurogenomic response to courtship was much greater in the wasp *V. canescens*. The differences may rely on the contrasted mating systems in these species. Despite the cost of reproduction, females mate multiple times in the majority of animal species in the wild including *D. melanogaster* [72] and *D. pseudoobscura* [73], most often with different males [74–76]. The benefits of polyandry for females include an adequate sperm supply [77], an increase in sperm competition [74] and a reduction in the cost of inbreeding [7]. In species with sl-CSD sex determination, mate choice is particularly determinant for female fitness, given the risk of genetic incompatibility. Furthermore, in monandrous species, such as *V. canescens* and a majority of parasitoid wasps (80%) [78], all progeny will have the same genitor. Thus, it is likely that selective pressures regarding mate choice should be stronger compared to those of polyandrous females. We suggest that the mating system might be an important determinant influencing the extent of neurogenomic response to courtship.

From a functional point of view, transcription factor activity and reproductive behaviour were some of the functions regulated following courtship. Indeed, we identified numerous transcription factors regulated following courtship, consistent with the large number of DEGs observed. In particular, we emphasized the transcription factors orthologous to *Thd1* and *muscleblind* and the kinase *Pink1*. All of these genes were also differentially expressed following courtship in *D. melanogaster* [45]. These candidates may have a conserved regulation pattern following courtship and might be associated with the response to courtship in insects. Our results also indicated the regulation of genes related to neurotransmitters, such as dopamine. The dopamine transporter *DAT* is overexpressed in courted females and the *ebony* gene involved in dopamine catabolism is downregulated in courted females, which is compatible with an increase in dopamine concentration following courtship. Dopamine is notably implicated in the control of motivation, movement and memory in the fruit fly *D. melanogaster* [79], and increased dopamine levels result in increased responsiveness to courtship cues [80]. Furthermore, we detected three neuropeptides downregulated in females courted by an unrelated male: *Capability*, *PDF* and *SIFamide*. These neuropeptides are involved in female sexual behaviour in the fruit fly, since reduction or absence of *SIFamide* makes females extremely receptive [81], while *PDF*-mutant females show an increased frequency of re-mating compared to wild-type females [82]. These three neuropeptides are candidates for involvement in the modulation of female receptivity. Together, these results demonstrate that the neurogenomic response to courtship in *V. canescens* involves neurotransmitters and neuropeptides. These genes are prime targets for further functional analyses.

Although courtship accounted for the majority of the detected DEGs, our experimental design nonetheless highlighted the major influence of relatedness between partners on female response to courtship. In total, 9.1% of the transcriptome was differentially expressed when comparing females courted by related males to those courted by unrelated males. In addition to the 2,780 courtship-responding genes that exhibited the same expression pattern whatever the relatedness of the courting male (modules 7 and 8), the network analysis highlighted 3 modules of sib-responsive genes (modules 4, 5 and 6; 2,824 genes) and 3 modules of genes regulated only in females courted by unrelated males (modules 1, 2 and 3; 1,239 genes) that could be related to female receptivity. Our results clearly showed that the relatedness of the courting male had a
major influence on the female response to courtship and that sib mating avoidance behaviour
observed in this species is correlated to complex and massive changes in gene expression.

Concerning the functions associated with genes that vary according to the relatedness of
the courting male, the GO terms Reproductive Behaviour and Male Mating Behaviour are
notable. All genes underlying both terms are homologous to the *yellow—major royal jelly proteins* (*yellow-mrjp*). The *yellow* gene is unique to insects and some bacteria, while the *mrjp* gene family is restricted to Hymenopteran species and evolved from recent duplications of the *yellow* gene [83]. In the honey bee, the *mrjp* gene family is involved in both group social behaviour (royal jelly is constituted with 90% of MRJP proteins), and in individual sexual behaviour, with sex- and caste-specific gene expression patterns [81]. *Yellow-mrjp* functions in parasitic wasps are unknown, even though the largest expansion of the family described so far came from the Nasonia parasitic wasp genome, where genes are expressed broadly in different tissues and life stages [84]. In this study, we identified 6 members of the *yellow-mrjp* family that were DEGs following courtship, of which one presented differential expression according to the relatedness of the courting male. Further functional characterization will be required to test whether the *yellow-mrjp* gene family is involved in sib mating avoidance and female receptivity. Among the regulated genes, we also highlighted transcription factors such as *PAX6* that could drive the transcriptomic changes accompanying female mate choice. Very few studies have explored the molecular pathways underlying kin recognition. In the amphibian *Xenopus laevis*, tadpoles exhibited plasticity in social preferences according to exposure to kinship odourants [85]. Sustained kin odourants exposure during development drives changes in neurotransmitter expression from GABA to dopamine neurons, which are stimulated from an increase in the expression of the transcription factor *PAX6* and accompanied by a behavioural preference for kin odourants [85]. Here, we observed that the *V. canescens* *PAX6* orthologue is downregulated in females courted by brothers. Further studies to characterize *PAX6* function in *V. canescens*, particularly in the context of kin recognition, could test whether common molecular pathways could be elicited for kin recognition in distant taxa such as amphibians and insects.

We had formulated two non-mutually exclusive hypotheses. First, the perception of courtship was mediated by a change in gene expression, that would result in similar expression patterns in all females regardless of their relationship to the courting male. We identified such patterns for 2,780 genes (modules 7 and 8). Second, changes in female receptivity could result in changes in transcriptomic profiles. In this case, similar expression patterns would be expected for isolated females and females courted by their brothers. We have identified such expression patterns for 1,239 genes (modules 1, 2 and 3). Thus our results suggest that both courtship perception and changes in female receptivity induce a different neurogenomic response. In addition to sib mating avoidance, *V. canescens* females express kin recognition in the context of host choice, since females prefer to lay their eggs in a host parasitized by others than by a relative [48]. Neurogenomic analysis of responses to the presence of a relative in distinct ecological contexts would determine whether there are molecular bases associated with kin recognition.

**Supporting information**

S1 Fig. Exclusion of one biological replicate from the analysis due to outlier status. Multi-variate analysis based on the expression profiles of the 500 genes with the highest variance across all samples showed that one biological replicate corresponding to females courted by related males (2 empty circles representing 2 technical replicates) is far from the other points (full circles) A) in the plan defined by the two first axes of the principal component analysis and B) in the sample hierarchical clustering dendrogram. Significant statistical support of
outlier status is indicated by an approximately unbiased *p*-value with one million multiscale bootstrap replicates.
(TIF)

**S2 Fig.** Boxplot showing expression patterns of the 6 members of the *yellow-mrjp* family were regulated following courting displays and/or according to the relatedness between the **female and the courting male**. Boxplot colours indicate biological condition: grey, isolated females; blue, females courted by unrelated males; orange, females courted by related males. The Y axes show the normalized counts after log transformation and VST normalization, and the boxplot whiskers show the range of reads between biological replicates. *, *P*-adj (FDR) < 0.01; NS, not significant.
(TIF)

**S3 Fig.** Boxplots showing the expression patterns of the top 5 hub genes for each of the 11 modules varied according to experimental conditions. The Y axes show the normalized counts after log transformation and VST normalization, and the boxplot whiskers show the range of reads between biological replicates.
(TIF)

**S1 Table.** Summary statistics of the sequencing runs and transcriptome assembly.
(XLSX)

**S2 Table.** Control females versus females courted by both related and unrelated males: A) list of 1,001 DEGs; B) 7 GO terms enriched in control females; C) 36 GO terms enriched in courted females.
(XLSX)

**S3 Table.** Females courted by related males versus females courted by unrelated males: A) list of 831 DEGs; B) 22 GO terms enriched in females courted by unrelated males; C) 66 GO terms enriched in females courted by related males.
(XLSX)

**S4 Table.** Enriched GO terms associated with the 11 significant modules of co-expressed genes.
(XLSX)

**S1 File.**
(CSV)

**S2 File.**
(CSV)

**S3 File.**
(TXT)

**S4 File.**
(CSV)

**S5 File.**
(R)

**Acknowledgments**

The genome sequences of *Venturia canescens* were kindly provided by Dr. JM Drezen (IRBI, Univ. Tours, France). This work was performed using the computing facilities of the CC LBBE/PRABI.
## Author Contributions

**Conceptualization:** Aurore Gallot, Emmanuel Desouhant.

**Data curation:** Aurore Gallot.

**Formal analysis:** Aurore Gallot.

**Funding acquisition:** Aurore Gallot, Emmanuel Desouhant.

**Investigation:** Aurore Gallot, Sandrine Sauzet.

**Methodology:** Emmanuel Desouhant.

**Project administration:** Emmanuel Desouhant.

**Visualization:** Aurore Gallot.

**Writing – original draft:** Aurore Gallot.

**Writing – review & editing:** Aurore Gallot, Emmanuel Desouhant.

## References

1. Andersson MB. Sexual selection. Princeton University Press; 1994.

2. Mays H Jr, Hill G. Choosing mates: good genes versus genes that are a good fit. *Trends Ecol Evol.* 2004; 19(10):554–9. https://doi.org/10.1016/j.tree.2004.07.018 PMID: 16701321

3. Ryan MJ, Cummings ME. Perceptual Biases and Mate Choice. Annu Rev Ecol Evol Syst. 2013; 44 (1):437–59.

4. Fuller RC, Houle D, Travis J. Sensory bias as an explanation for the evolution of mate preferences. Am Nat. 2005; 166(4):437–46. https://doi.org/10.1086/444443 PMID: 16224700

5. Møller AP, Alatalo RV. Good-genes effects in sexual selection. Proc R Soc Lond B Biol Sci. 1999; 266 (1414):85.

6. Neff BD, Pitcher TE. Genetic quality and sexual selection: an integrated framework for good genes and compatible genes. Mol Ecol. 2005; 14(1):19–38. https://doi.org/10.1111/j.1365-294X.2004.02395.x PMID: 15643948

7. Tregenza T, Wedell N. Genetic compatibility, mate choice and patterns of parentage: invited review. Mol Ecol. 2000; 9(8):1013–27. https://doi.org/10.1046/j.1365-294x.2000.00964.x PMID: 10964221

8. Pusey A, Wolf M. Inbreeding avoidance in animals. *Trends Ecol Evol.* 1996; 11(5):201–6. https://doi.org/10.1016/0169-5347(96)0028-8 PMID: 21237809

9. Rosengarten RD, Nicotra ML. Model Systems of Invertebrate Allorecognition. Curr Biol. 2011; 21(2): R82–92. https://doi.org/10.1016/j.cub.2010.11.061 PMID: 21256442

10. Taketa DA, De Tomaso AW. Botryllus schlosseri allorecognition: tackling the enigma. Dev Comp Immunol. 2015; 48(1):254–65. https://doi.org/10.1016/j.dci.2014.03.014 PMID: 24709050

11. Beauchamp GK, Yamazaki K, Boyse EA. The chemosensory recognition of genetic individuality. Sci Am. 1985; 253(1):86–93. https://doi.org/10.1038/scientificamerican0785-86 PMID: 4001916

12. Boyse EA, Beauchamp GK, Yamazaki K. The genetics of body scent. *Trends Genet.* 1987; 3:97–102.

13. Potts WK, Manning CJ, Wakeland EK. Mating patterns in seminatural populations of mice influenced by MHC genotype. *Nature*. 1991; 352(6336):619. https://doi.org/10.1038/352619a0 PMID: 1865924

14. Singh PB, Brown RE, Roser B. MHC antigens in urine as olfactory recognition cues. *Nature*. 1987; 327 (6118):161–4. https://doi.org/10.1038/327161a0 PMID: 3574476

15. Ben-Shaul Y, Katz LC, Mooney R, Dulac C. In vivo vomeronasal stimulation reveals sensory encoding of conspecific and allospecific cues by the mouse accessory olfactory bulb. Proc Natl Acad Sci U S A. 2010; 107(11):5172–7. https://doi.org/10.1073/pnas.0915147107 PMID: 20194746

16. Chamero P, Leinders-Zufall T, Zufall F. From genes to social communication: molecular sensing by the vomeronasal organ. *Trends Neurosci.* 2012; 35(10):597–606. https://doi.org/10.1016/j.tins.2012.04.011 PMID: 22658923

17. Manning C, Wakeland E, Potts W. Communal nesting patterns in mice implicate MHC genes in kin recognition. *Nature*. 1992; 360(6404):581–3. https://doi.org/10.1038/360581a0 PMID: 1461279
18. Milinski M, Griffiths S, Wegner KM, Reusch TBH, Haas-Assenbaum A, Boehm T. Mate choice decisions of stickleback females predictably modified by MHC peptide ligands. Proc Natl Acad Sci U S A. 2005; 102(12):4414. https://doi.org/10.1073/pnas.0408264102 PMID: 15755811

19. Yamazaki K, Beauchamp GK, Singer A, Bard J, Boyse EA. Odortypes: Their origin and composition. Proc Natl Acad Sci. 1999; 96(4):1522. https://doi.org/10.1073/pnas.96.4.1522 PMID: 9990056

20. Derstine NT, Gries R, Zhai H, Jimenez SI, Gries G. Cuticular hydrocarbons determine sex, caste, and nest membership in each of four species of yellowjackets (Hymenoptera: Vespidae). Insectes Sociaux. 2018; 65(4):581-91.

21. Howard RW, Blomquist GJ. Ecological, behavioral, and biochemical aspects of insect hydrocarbons. Annu Rev Entomol. 2005; 50:371-93. https://doi.org/10.1146/annurev.ento.50.071803.130359 PMID: 15355247

22. Singer TL. Roles of hydrocarbons in the recognition systems of insects. Am Zool. 1998; 38(2):394-405.

23. Van Wilgenburg E, Symonds MRE, Elgar MA. Evolution of cuticular hydrocarbon diversity in ants. J Evol Biol. 2011; 24(6):1188-98. https://doi.org/10.1111/j.1420-9101.2011.02248.x PMID: 21375648

24. Couto A, Mitra A, Thiery D, Marion-Poll F, Sandoz J-C. Hornets Have It: A Conserved Olfactory Subsystem for Social Recognition in Hymenoptera? Front Neuroanat. 2017; 11:48. https://doi.org/10.3389/fnana.2017.00048 PMID: 28659767

25. Kelber C, Roessler W, Kleineidam CJ. Phenotypic plasticity in number of glomeruli and sensory innervation of the antennal lobe in leaf-cutting ant workers (A. vollenweideri). Dev Neurobiol. 2010; 70 (4):222-34. https://doi.org/10.1002/dneu.20782 PMID: 20029932

26. McKenzie SK, Fetter-Pruneda I, Ruta V, Kronauer DJC. Transcriptomics and neuroanatomy of the clonal raider ant implicate an expanded clade of odorant receptors in chemical communication. Proc Natl Acad Sci U S A. 2016; 113(49):14091-6. https://doi.org/10.1073/pnas.1610800113 PMID: 27911792

27. Nishikawa M, Watanabe H, Yokohari F. Higher brain centers for social tasks in worker ants, Camponotus japonicus. J Comp Neurol. 2012; 520(7):1584-98. https://doi.org/10.1002/cne.23001 PMID: 22102363

28. Zube C, Roessler W. Caste- and sex-specific adaptations within the olfactory pathway in the brain of the ant Camponotus floridanus. Arthropod Struct Dev. 2008; 37(6):469-79. https://doi.org/10.1016/j.asd.2008.05.004 PMID: 18621145

29. Cummings ME. The mate choice mind: studying mate preference, aversion and social cognition in the female poeciliid brain. Anim Behav. 2015; 103:249-58.

30. Bloch NI, Corral-López A, Buechel SD, Kotrschal A, Kolm N, Mank JE. Early neurogenomic response associated with variation in guppy female mate preference. Nat Ecol Evol. 2018; 2(11):1772-81. https://doi.org/10.1038/s41559-018-0682-4 PMID: 30297748

31. Cummings ME, Larkins-Ford J, Reilly CRL, Wong RY, Ramsey M, Hofmann HA. Sexual and social stimuli elicit rapid and contrasting genomic responses. Proc R Soc B-Biol Sci. 2008; 275(1633):393-402.

32. Ramsey ME, Maginnis TL, Wong RY, Brock C, Cummings ME. Identifying context-specific gene profiles of social, reproductive, and mate preference behavior in a fish species with female mate choice. Front Neurosci. 2012; 6:UNSP 62.

33. van Wilgenburg E, Driessen G, Beukeboom LW. Single locus complementary sex determination in Hymenoptera: an « unintelligent » design? Front Zool. 2006; 3:1. https://doi.org/10.1186/1742-9994-3-1 PMID: 16393347

34. Cook JM. Sex determination in the Hymenoptera: a review of models and evidence. Heredity. 1993; 71:421.

35. Ross KG, Vargo EL, Keller L, Trager JC. Effect of a founder event on variation in the genetic sex-determining system of the fire ant Solenopsis invicta. Genetics. 1993; 135(3):843-54. PMID: 8293983

36. Heimpel GE, de Boer JG. Sex determination in the hymenoptera. Annu Rev Entomol. 2008; 53:209-30. https://doi.org/10.1146/annurev.ento.53.103106.093441 PMID: 17803453

37. Cowan DP, Stahlhut JK. Functionally reproductive diploid and haploid males in an inbreeding hymenopteran with complementary sex determination. Proc Natl Acad Sci U S A. 2004; 101(28):10374-9. https://doi.org/10.1073/pnas.0402481101 PMID: 15232002

38. Elias J, Mazzi D, Dorn S. No Need to Discriminate? Reproductive Diploid Males in a Parasitoid with Complementary Sex Determination. PLOS ONE. 2009; 4(6):e6024. https://doi.org/10.1371/journal.pone.0006024 PMID: 19551142

39. Zayed A, Packer L. Complementary sex determination substantially increases extinction proneness of haplodiploid populations. Proc Natl Acad Sci U S A. 2005; 102(30):10742-6. https://doi.org/10.1073/pnas.0502271102 PMID: 16020532
40. Beukeboom LW. Single-locus complementary sex determination in the ichneumonid *Venturia canescens* (Gravenhorst)(Hymenoptera). Neth J Zool. 2001; 51(1):1-15.

41. Metzger M, Fischbein D, Auguste A, Fauvergue X, Bernstein C, Desouhant E. Synergy in information use for mate finding in a parasitic wasp. J Anim Ecol. 2010;

42. Metzger M, Bernstein C, Hoffmeister TS, Desouhant E. Does Kin Recognition and Sib-Mating Avoidance Limit the Risk of Genetic Incompatibility in a Parasitic Wasp? Plos One. 2010; 5(10):e13505.
https://doi.org/10.1371/journal.pone.0013505 PMID: 20976063

43. Zayed A, Robinson GE. Understanding the Relationship Between Brain Gene Expression and Social Behavior: Lessons from the Honey Bee. In: Bassler BL, editor. Annual Review of Genetics, Vol 46. 2012. p. 591-615. https://doi.org/10.1146/annurev-genet-110711-155517 PMID: 22994354

44. Van Santen M, Schneider M. Courtship description in the solitary parasitoid wasp *Venturia canescens* (Hymenoptera). Proc Exp Appl Entomol. 2002; 13:87-92.

45. Immonen E, Ritchie MG. The genomic response to courtship song stimulation in female *Drosophila melanogaster*. Proc R Soc B-Biol Sci. 2012; 279(1732):1359-65.

46. Veltsos P, Fang Y, Cossins AR, Snook RR, Ritchie MG. Mating system manipulation and the evolution of sex-biased gene expression in Drosophila. Nat Commun. 2017; 8(1):2072. https://doi.org/10.1038/s41467-017-02232-6 PMID: 29233985

47. Salt G. The hosts of *Nemeritis canescens* a problem in the host specificity of insect parasitoids. Ecol Entomol. 1976; 1(1):63-7.

48. Amat I, Desouhant E, Bernstein C. Differential use of conspecific-derived information by sexual and asexual parasitic wasps exploiting partially depleted host patches. Behav Ecol Sociobiol. 2009; 63 (4):563-72.

49. Chuine A, Sauzet S, Debias F, Desouhant E. Consequences of genetic incompatibility on fitness and mate choice: the male point of view. Biol J Linn Soc. 2015; 114(2):279-86.

50. Rittschof CC, Robinson GE. Chapter Five—Behavioral Genetic Toolkits: Toward the Evolutionary Origins of Complex Phenotypes. In: Orgogozo V, editor. Current Topics in Developmental Biology. Academic Press; 2016. p. 157-204.

51. Todd EV, Black MA, Gemmell NJ. The power and promise of RNA-seq in ecology and evolution. Mol Ecol. 2016; 25(6):1224-41. https://doi.org/10.1111/mec.13526 PMID: 26756714

52. Leobold M, Bézier A, Pichon A, Herniou EA, Volkoff A-N, Drezen J-M. The Domestication of a Large DNA Virus by the Wasp *Venturia canescens* Involves Targeted Genome Reduction through Pseudo-genization. Genome Biol Evol. 1 juil 2018; 10(7):1745-64. https://doi.org/10.1093/gbe/evy127 PMID: 29931159

53. Pichon A, Bézier A, Urbach S, Aury J-M, Jouan V, Ravallec M, et al. Recurrent DNA virus domestication leading to different parasite virulence strategies. Sci Adv. 2015; 1(10).

54. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, et al. Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks. Nat Protoc. 2012; 7(3):562-78.
https://doi.org/10.1038/nprot.2012.016 PMID: 22383036

55. Anders S, Pyl PT, Huber W. HTSeq a Python framework to work with high-throughput sequencing data. Bioinforma Oxf Engl. 2015; 31(2):166-9.

56. R Core Team. R: A language and environment for statistical computing [Internet]. R Fondation for Statistical Computing, Austria; 2017.

57. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014; 15(12):550. https://doi.org/10.1186/s13059-014-0550-8 PMID: 25516281

58. Suzuki R, Shimodaira H. Pvclust: an R package for assessing the uncertainty in hierarchical clustering. Bioinformatics. 2006; 22(12):1540-2. https://doi.org/10.1093/bioinformatics/btl117 PMID: 16595560

59. Benjamini Y, Hochberg Y. Controlling the False Discovery Rate—a practical and powerful approach. J R Stat Soc Ser B-Methodol. 1995; 57(1):289-300.

60. Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008; 9:559. https://doi.org/10.1186/1471-2105-9-559 PMID: 19114008

61. Langfelder P, Mischel PS, Horvath S. When Is Hub Gene Selection Better than Standard Meta-Analysis? PLOS ONE. 2013; 8(4):e61505. https://doi.org/10.1371/journal.pone.0061505 PMID: 23613865

62. Filteau M, Pavey SA, St-Cyr J, Bernatchez L. Gene Coexpression Networks Reveal Key Drivers of Phenotypic Divergence in Lake Whitefish. Mol Biol Evol. 2013; 30(6):1384-96. https://doi.org/10.1093/molbev/mst053 PMID: 23519315
63. Conesa A, Gotz S, Garcia-Gomez JM, Terol J, Talon M, Robles M. Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005; 21(18):3674-6. https://doi.org/10.1093/bioinformatics/bti610 PMID: 16081474

64. McNeill MS, Kapheim KM, Brockmann A, McGill TAW, Robinson GE. Brain regions and molecular pathways responding to food reward type and value in honey bees. Genes Brain Behav. 2016; 15(3):305-17. https://doi.org/10.1111/gbb.12275 PMID: 26566901

65. Juni N, Yamamoto D. Genetic analysis of chaste, a new mutation of *Drosophila melanogaster* characterized by extremely low female sexual receptivity. J Neurogenet. 2009; 23(3):329-40. https://doi.org/10.1080/01677060802471601 PMID: 19169922

66. Ferguson LC, Green J, Surridge A, Jiggins CD. Evolution of the Insect Yellow Gene Family. Mol Biol Evol. 2011; 28(1):257-72. https://doi.org/10.1093/molbev/msq192 PMID: 20656794

67. Callaerts P, Leng S, Clements J, Benassayag C, Cribbs D, Kang YY, et al. Drosophila *Pax-6/eyeless* is essential for normal adult brain structure and function. J Neurobiol. 2001; 46(2):73-88. https://doi.org/10.1002/1097-4695(20010205)46:2<73::aid-neu10>3.0.co;2-n PMID: 11153010

68. Chandrasekaran S, Ament SA, Eddy JA, Rodriguez-Zas SL, Schatz BR, Price ND, et al. Behavior-specific changes in transcriptional modules lead to distinct and predictable neurogenomic states. Proc Natl Acad Sci U S A. 2011; 108(44):18020-5. https://doi.org/10.1073/pnas.1114093108 PMID: 21960440

69. Whitney O, Pfenning AR, Howard JT, Blatti CA, Liu F, Ward JM, et al. Core and region-enriched networks of behaviorally regulated genes and the singing genome. Science. 2014; 346(6215):1334.

70. Clayton DF. The genomic action potential. Neurobiol Learn Mem. 2000; 74(3):185-216. https://doi.org/10.1006/nlme.2000.3967 PMID: 11031127

71. Saha RN, Wissink EM, Bailey ER, Zhao M, Fargo DC, Hwang J-Y, et al. Rapid activity-induced transcription of Arc and other IEGs relies on poised RNA polymerase II. Nat Neurosci. 2011; 14(7):848-U62. https://doi.org/10.1038/nn.2839 PMID: 21623364

72. Imhof M, Harr B, Brem G, Schlötterer C. Multiple mating in wild *Drosophila melanogaster* revisited by microsatellite analysis. Mol Ecol. 1998; 7(7):915-7. https://doi.org/10.1046/j.1365-294x.1998.00382.x PMID: 9691492

73. Anderson WW. Frequent multiple insemination in a natural population of Drosophila pseudoobscura. Am Nat. 1974; 108(963):709-11.

74. Simmons L. The evolution of polyandry: patterns of genotypic variation in female mating frequency, male fertilization success and a test of the sexy-sperm hypothesis. J Evol Biol. 2003; 16(4):624-34. https://doi.org/10.1046/j.1420-9101.2003.00572.x PMID: 14632226

75. Hosken D, Stockley P. Benefits of polyandry: a life history perspective. In: Evolutionary Biology. Springer; 2003. p. 173-94.

76. Simmons LW. The evolution of polyandry: sperm competition, sperm selection, and offspring viability. Annu Rev Ecol Evol Syst. 2005; 36:125-46.

77. Turner ME, Anderson WW. Multiple mating and female fitness in *Drosophila pseudoobscura*. Evol Int J Org Evol. 1983; 37(4):714-23.

78. Ridley M. Clutch size and mating frequency in parasitic Hymenoptera. Am Nat. 1993; 142(5):893-910.

79. Waddell S. Dopamine reveals neural circuit mechanisms of fly memory. Trends Neurosci. 2010; 33 (10):457-64. https://doi.org/10.1016/j.tins.2010.07.001 PMID: 20701984

80. Andretic R, van Swinderen B, Greenspan RJ. Dopaminergic modulation of arousal in Drosophila. Curr Biol. 2005; 15(13):1165-75. https://doi.org/10.1016/j.cub.2005.05.025 PMID: 16005288

81. Terhzaz S, Rosay P, Goodwin SF, Veenstra JA. The neuropeptide SIFamide modulates sexual behavior in Drosophila. Biochem Biophys Res Commun. 2007; 352(2):305-10. https://doi.org/10.1016/j.bbrc.2006.11.030 PMID: 17126293

82. Krupp JJ, Billeter J-C, Wong A, Choi C, Nitabach MN, Levine JD. Pigment-dispersing factor modulates pheromone production in clock cells that influence mating in drosophila. Neuron. 2013; 79(1):54-68. https://doi.org/10.1016/j.neuron.2013.05.019 PMID: 23849197

83. Drapeau MD, Albert S, Kucharski R, Prusko C, Maleszka R. Evolution of the Yellow/Major Royal Jelly Protein family and the emergence of social behavior in honey bees. Genome Res. 2006; 16(11):1385-94. https://doi.org/10.1101/gr.5012006 PMID: 17065613

84. Werren JH, Richards S, Desjardins CA, Niehuis O, Gadau J, Colbourne JK, et al. Functional and evolutionary insights from the genomes of three parasitoid Nasonia species. Science. 2010; 327(5963):343-8. https://doi.org/10.1126/science.1178028 PMID: 20075255

85. Dulcis D, Lippi G, Stark CJ, Do LH, Berg DK, Spitzer NC. Neurotransmitter Switching Regulated by miRNAs Controls Changes in Social Preference. Neuron. 13 2017; 95(6):1319. https://doi.org/10.1016/j.neuron.2017.08.023 PMID: 28867550