The molecular basis of socially-mediated phenotypic plasticity in a eusocial paper wasp

The molecular basis of socially-mediated phenotypic plasticity in a eusocial paper wasp | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article The molecular basis of socially-mediated phenotypic plasticity in a eusocial paper wasp Benjamin Taylor, Alessandro Cini, Christopher Wyatt, Max Reuter, and 1 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-43815/v1 This work is licensed under a CC BY 4.0 License Status: Published Journal Publication published 02 Feb, 2021 Read the published version in Nature Communications → Version 1 posted You are reading this latest preprint version Abstract Phenotypic plasticity, the ability to produce multiple phenotypes from a single genotype, represents an excellent model with which to examine the relationship between gene expression and phenotypes. Despite this, analyses of the molecular bases of plasticity have been limited by the challenges of linking individual phenotypes with individual-level gene expression profiles, especially in the case of complex social phenotypes. Here, we tackle this challenge by analysing the individual-level gene expression profiles of Polistes dominula paper wasps following the loss of a queen, a perturbation that induces some individuals to undergo a significant phenotypic shift and become replacement reproductives. Using a machine learning approach, we find a strong response of caste-associated gene expression to queen loss, wherein individuals’ expression profiles become intermediate between queen and worker states. Importantly, this change occurs even in individuals that appear phenotypically unaffected. Part of this response is explained by individual attributes, most prominently age. These results demonstrate that large changes in gene expression may occur in the absence of detectable phenotypic changes, resulting here in a socially mediated de-differentiation of individuals at the transcriptomic but not the phenotypic level. Our findings also highlight the complexity of the relationship between gene expression and phenotype, where transcriptomes are neither a direct reflection of the genotype nor a proxy for the molecular underpinnings of the external phenotype. Entomology Animal Behavior Bioinformatics Phenotypic plasticity machine learning Polistes dominula Figures Figure 1 Figure 2 Figure 3 Introduction The relationship between gene expression and external phenotype is complex and unresolved. Much research in behavioural and evolutionary ecology is based on the implicit assumption that phenotypic traits can be modelled as though they directly reflect gene expression patterns, and that evolutionary trajectories can therefore be studied while remaining agnostic with regard to the underlying molecular mechanisms (Hadfield et al 2007 ; Rittschof & Robinson 2014 ; Rubin 2016 ). This ‘phenotypic gambit’ has proven a useful rule of thumb, permitting the establishment of a rich body of literature surrounding the evolution of complex traits, despite a lack of data relating to the genetic basis of these traits (e.g. Réale et al 2010 ; Chapman et al 2011 ; Fowler-Finn & Rodríguez 2012 ). In the past decade, however, advances in the affordability of ‘omic’ data and availability of powerful bioinformatic methods have greatly enhanced our ability to assess the assumptions made by the phenotypic gambit (Rittschof & Robinson 2014 ; Heyes 2016 ). The time is right to disentangle the molecular bases of complex phenotypic traits. Phenotypic plasticity, the ability of an individual to effect phenotypic changes in response to external cues, is an ideal phenomenon with which to study the relationship between gene expression and phenotype because it involves the production of multiple phenotypes without gene sequence changes. Of particular value are species in which adult individuals can be experimentally induced to transition between distinct, measurable phenotypes. By comparing the gene expression profiles of groups of individuals that differ in their phenotypes as a result of plasticity, it is possible to isolate phenotypic effects of gene expression while controlling for genetic differences. Using this approach, significant progress has been made in unravelling the molecular underpinnings of sequential sex changes in hermaphroditic fish (Horiguchi et al 2013 ; Casas et al 2016 ; Todd et al 2019 ), the distinct gregarious social phenotype of desert locusts (Cullen et al 2017 ; Lo et al 2018 ), and the reproductive castes of social insects (Simola et al 2016 ; Libbrecht et al 2018 ; Rehan et al 2018 ; Shell & Rehan 2018 ). Such studies typically rely on comparisons between groups of individuals with well-differentiated phenotypes, however. As a result, little is known about more subtle effects during the transition from one morph to the other, and the relationship between expression patterns and phenotypic traits at the individual level. The reproductive castes found in the colonies of social insects provide excellent model systems for determining the extent to which fine-scale changes in phenotype are reflected at the molecular level. With a few exceptions (reviewed in Schwander et al 2010 ), the distinct queen and worker phenotypes found in such colonies are plastically determined either during development or in adulthood. Workers in some species can be experimentally induced to transition to a reproductive role in response to the removal of a colony’s queen (e.g. Strassmann et al 2004 ; Tibbetts & Huang 2010 ) or as a result of exposure to varying levels of brood (e.g. Chandra et al 2018 ; Libbrecht et al 2018 ), allowing changes in the behavioural, physiological and molecular traits that define caste identity to be tracked. An additional benefit—and challenge—of studying social insect colonies is that they involve complex social structures. Such interactions can be hard to study, but offer the opportunity to assess the effects of social interactions upon phenotypes and transcriptomes. In the present study, we explore the relationship between transcriptome and phenotype in the context of plastic caste expression in the European paper wasp Polistes dominula (Christ 1791), which is a model organism often used in studies of social insect behaviour (reviewed in Starks & Turillazzi 2006; Jandt et al 2014 ) and, more recently, for analyses of caste gene expression (e.g. Patalano et al 2015 ; Standage et al 2016 ; Geffre et al 2017 ; Manfredini et al 2018 ). In this species, removing the established queen from a single-foundress colony induces a queen succession process in which one (or very few) workers transition to a queen phenotype, with age playing a key role in predicting which individual will do so (Pardi 1948 ; Tsuji &Tsuji 2005 ): almost invariably, the new queen is one of the oldest individuals, and there is little conflict over succession (Strassmann et al 2004 ; Taylor et al 2020 ). In a recent paper (Taylor et al 2020 ) we followed responses to queen removal in P. dominula on a fine scale by measuring individual-level behavioural and physiological traits and generating a univariate measure of caste identity (‘queenness’) to describe individuals’ phenotypic profiles ( Box 1 ). We found that phenotypic responses were limited to the subset of individuals that transitioned to the queen role, while other individuals exhibited little or no measurable change in measured traits. This system, in which individuals within a controlled environment vary strongly and predictably in their phenotypic response to a shared stimulus, affords an excellent opportunity to track the relationship between gene expression and phenotypic expression at the individual level. To do so, here we complement the individual-level phenotypic data of Taylor et al. ( 2020 ) by assaying the brain gene expression profiles of 107 individuals for which fine-scale behavioural and ovarian data have been obtained, including queens and workers from stable colonies and individuals from colonies that had their queens experimentally removed. We then assess the relationship between the transcriptome and the phenotypic changes induced by queen removal. To match the univariate measure of caste identity that we applied to our phenotypic data ( Box 1 ), we use a support vector classification approach to generate a univariate measure of caste-associated gene expression. Support vector classification is a powerful tool with which to transform complex patterns in multidimensional data into a continuous classification score, allowing the detection of subtle, widespread signals of differential expression between phenotypic states that are likely to be missed in conventional differential expression analyses. Support vector machines (SVMs) have become a key tool in the early identification of phenotypically indistinguishable cancer subtypes (Segal et al 2003 ; Abeel et al 2010 ; Huang et al 2018 ) and their value has recently been demonstrated in animal behaviour studies (e.g. Chakravarty et al 2019 ; Rast et al 2020 ). Given their sensitivity, SVMs should be ideally suited to quantify gene expression variation across the spectrum between differentiated worker and queen roles, but to our knowledge this is the first time SVMs have been used to interrogate the molecular basis of behavioural traits. Applying an SVM approach along with standard differential expression and gene co-expression analyses, we show that brain gene expression responses to queen removal in P. dominula include a colony-wide response that does not match that observed at the phenotypic level. Our results indicate that gene expression in P. dominula colonies reflects both a generalised response to queen loss which is independent of phenotype, and a phenotype-specific response that tracks individuals’ expression of plastic phenotypic changes. To our knowledge, this study provides the most comprehensive analysis to date of the ways that plastic phenotypes are reflected at the transcriptomic level; our results expose the complexity of the relationship between individual-level gene expression and individual phenotypes and call into question the assumptions of the phenotypic gambit. Results Support vector classification reveals consistent patterns of caste gene expression differentiation involving many genes Our support vector approach allowed us to reduce complex variation in brain gene expression data down to a single dimension of caste identity, matching the Bayesian logistic regression model we used in Taylor et al. ( 2020 ) to condense ovarian and behavioural data into a unidimensional metric of phenotypic caste (‘queenness’; Box 1 ). We trained an SVM using gene expression data from 27 queens and 12 workers from stable, unmanipulated colonies. An initial model, based on all 10734 genes annotated in our experiment, achieved a root mean squared validation error of 0.065 in three-fold cross-validation, i.e. a model trained on a random subset of two thirds of workers and queens classified the remaining third within 25.5% of their true values (worker = 0, queen = 1). By applying a process of iterative feature selection we were able to remove genes that had low weights in the classification and contributed to overfitting ( Figure S1A; Table S1 ). Doing so allowed us to identify an optimal model containing 1992 genes with a substantially reduced root mean squared classification error of 0.021 ( Table S2 ). This model classified queens and workers very consistently, with strong separation of queens from workers (Fig. 1A). Thirty-four gene ontology (GO) terms were significantly enriched among these 1992 genes, including a number of terms associated with translation such as rRNA processing , tRNA aminoacylation , and ribosomal large subunit biogenesis ( Figure S1B ). The 1992 genes in our optimised SVM represent a much larger set than that which is found to differentiate queens and workers based on a standard differential expression approach. When we applied a DESeq2 analysis with a 1.5 fold-change threshold to the same set of queens and workers used to train the SVM, we identified just 81 differentially expressed genes (with no associated GO terms), a small number that is typical for similar analyses in Polistes (e.g. Toth et al 2014 ; Patalano et al 2015 ; Geffre et al 2017 ). All but four of these genes were present in the larger SVM set ( Figure S1C ; Table S3 ), suggesting that the SVM captures almost all the information obtained from these standard methods. Conversely, an SVM trained using just these 81 differentially-expressed genes exhibited a cross-validation error rate of 0.065, no better than the original un-optimised model. Thus, the picture of caste differentiation provided by standard differential expression analysis appears to miss a great number of subtle differences in individuals’ gene expression profiles that contribute to caste differentiation. Colony-wide brain gene expression responses to queen removal Following the loss of a queen from a Polistes colony, typically one or a few individuals undergo a phenotypic transition to become a replacement queen while the rest of the colony members remain workers (Miyano 1991 ; Dapporto et al 2005 ; Taylor et al 2020 ). If this response is reflected at the level of gene expression, then we would expect to observe significant changes in the transcriptional profiles of individuals that exhibit a phenotypic transition, but not in those that do not. Our individual-level expression data allowed us to track responses in 62 individuals’ gene expression profiles following queen removal across 19 nests. We analysed individuals’ gene expression profiles at three days after queen removal, when queen replacement is ongoing (n = 24), and at twelve days after queen removal, when succession is largely settled at the phenotypic level (Strassmann et al 2004 ; Taylor et al 2020 ; n = 38). Individuals for sequencing were selected to cover a wide range of phenotypes, including those that remained entirely worker-like, those that had transitioned to highly queen-like phenotypes, and those with intermediate phenotypes at the time of sampling ( Figure S1D-F ). Analysing shifts in SVM estimates of individual wasps allowed us to assess the degree to which these concurred with the changes visible at the phenotypic level. Doing so, we found that the SVM estimates of individuals from post-removal nests were intermediate between those of queens and workers from control nests (Fig. 1B), a finding which concurs with the placement of these individuals according to principal component analysis ( Figure S1G ). This result is surprising given that the majority of individuals on queen removal nests are phenotypically indistinguishable from workers on control nests in terms of ovarian development and behavioural dominance (Taylor et al 2020 ). Thus, the large majority of individuals on queen removal nests exhibited perturbation of their caste-associated gene expression, even though only a few of these individuals exhibited responses to queen loss at the level of physiology or behaviour and transitioned from a worker to a queen role. Furthermore, we found no evidence that the degree of transcriptional perturbation declined over time following queen removal. SVM classification estimates of individuals from queen removal nests did not differ significantly between days three and twelve following queen removal (QR3 mean 0.331 ± 0.111; QR12 mean 0.345 ± 0.082; Wilcoxon W = 369, p = 0.55). The effects of queen removal on individuals’ gene expression profiles therefore appear to be both widespread and persistent, affecting all individuals in a nest and lasting beyond the point at which a new queen has already become phenotypically established. Interestingly, the strong colony-wide perturbation following queen loss that this SVM approach identifies would have been entirely missed using a standard differential expression approach: DESeq2 identified just five genes as differentially expressed between control workers and individuals from manipulated nests, with no associated GO enrichment ( Table S4 ). Age and queenness explain variation in individual-level molecular responses to queen loss While SVM classification indicates that queen removal causes colony-wide perturbation to brain expression profiles, classifier estimates varied substantially between individuals following queen removal, spanning a much greater range of values (0.116–0.540) than those of queens (0.900–1.100) or workers (0.057-0.100) from queenright colonies (Fig. 1). To better understand this variation, we examined whether the classifier estimates for individuals from manipulated colonies were predicted by those individuals' phenotypic traits—specifically ovarian development, behavioural dominance, and age ( Box 1 ; Taylor et al 2020 ). Our phenotypic measure of caste identity (queenness) was a significant predictor of expression-based SVM classifier estimates when fitted using a linear model (slope ± SE = 0.0414 ± 0.0114, p = 6.0 × 10 − 4 ; Fig. 2 A). The individual components of queenness, ovarian development and dominance were also significant or near-significant predictors of SVM classification individually (ovarian development: slope ± SE = 0.0426 ± 0.0113, p = 4.0 × 10 − 4 ; dominance: slope ± SE = 0.0231 ± 0.0122, p = 0.06). Notably, however, because phenotypic queenness was strongly correlated with age among post-removal individuals (cor = 0.4832, p = 8.0 × 10 − 5 ; Taylor et al 2020 ), age was an equally strong predictor of caste estimates when fitted in a separate linear model (slope ± SE = 0.0516 ± 0.0106, p = 9.8 × 10 − 6 ; Fig. 2 B). Thus, the significance of queenness as a predictor of caste estimates might have been an artefact of the fact that both are correlated with age. To test whether queenness had an effect over and above that accounted for by age, we calculated the residuals of phenotypic queenness on age and fitted the caste estimates against these. These residuals were not significantly predictive of SVM classification (slope ± SE = 0.0202 ± 0.0123, p = 0.11; Fig. 2 C), nor were the residuals of ovarian development alone on age (slope ± SE = 0.0223 ± 0.0123, p = 0.07) or the residuals of dominance alone (slope ± SE = 0.0114 ± 0.0126, p = 0.37). Age is thus the strongest determinant both of individuals’ caste phenotypes and of their caste-associated gene expression. In order to further discern the factors that shaped individuals’ gene expression profiles following queen removal, we performed weighted gene co-expression network analysis (WGCNA) using the full set of annotated genes. We generated 22 consensus modules across all samples, and then determined which traits significantly predicted the expression of a given module within each group (Fig. 2 D-F). In doing so, we identified three modules of particular interest. Module 11 ( Red ) consists of 519 genes associated with 14 enriched GO terms, including a number of terms associated with molecular binding (Fig. 2 G). This module is notable in that its expression in individuals from queen removal colonies is positively associated with SVM classifications but not with phenotypic correlates of caste. A second module, Module 2 ( Lightyellow ), consists of 102 genes and appears to be associated with many caste-related traits in the queen removal condition, being negatively correlated with queenness, age, SVM classification and (more weakly) with ovarian development among queen removal individuals. We also found weak evidence that this module is negatively associated with age among workers from control nests. Unexpectedly, despite the seeming importance of Module 2 in predicting caste identity among workers, not a single GO term was enriched among the genes in the module at p < 0.01. Finally, Module 8 ( Yellow ), consists of 614 genes and has a strongly negative correlation with age among both control workers and queen removal individuals. This module is associated with 21 GO terms, including a number of terms associated with molecular binding, such as protein binding, DNA binding, RNA binding and chromatin binding (Fig. 2 H). The latter term is particularly notable: regulation of chromatin accessibility is one of several epigenetic processes that have been implicated in the control of caste expression in social insects (Simola et al 2016 ; Wojciechowski et al 2018; Duncan et al 2020). Downregulation of this module with age might plausibly contribute to the increased phenotypic plasticity exhibited by older individuals in P. dominula. Discussion This study provides one of the most comprehensive analyses to date regarding the transcriptomic signatures of plastic phenotypes, and the first to employ an SVM approach to interrogate the molecular basis of behaviour. We have applied this approach to identify caste-specific gene expression profiles in P. dominula and to explore the relationship between transcriptomic and plastic phenotypic changes following a major social perturbation in paper wasp colonies. We identify a set of nearly 2000 genes that optimally capture gene expression differences between established queens and workers. Using a caste classifier based on these genes, we find that queen removal leads to a colony-wide shift in expression, where the expression profiles of all individuals move towards a state intermediate between those of established queens and workers. Individual variation around this intermediate state is related to age and phenotypic attributes, with older (and therefore more queen-like) individuals showing expression profiles that are closer to that expected for established queens. Our results show that molecular responses to queen removal in P. dominula consist of both a general colony-wide response independent of phenotypic change and a response that reflects the plastic phenotypic transition. Our study contributes to the enormous progress in our understanding of the relationship between molecular changes and changes in phenotypic expression that has been made in the past decade, facilitated by the increased availability of 'omic’ data and complex bioinformatic analyses. Recent studies have started to challenge the view that there is a direct correspondence between transcriptomic states and external phenotypes. Libbrecht et al ( 2018 ), for example, show that gene expression responses associated with a reversible phenotypic change differ qualitatively based on the directionality of the change (from reproductive to non-reproductive or vice versa). Meanwhile, molecular manipulations have revealed a surprising degree of plasticity in canonically implastic traits such as mammalian sex (Matson et al 2011 ) or ant castes (Simola et al 2016 ). Our results go further, showing a shift in caste-specific brain gene expression profiles among individuals whose phenotypic caste expression remains otherwise apparently unchanged. This shows that the expectation of a close match between expression profiles and phenotypes is excessively simplistic, or at least that detecting such a match requires detailed knowledge of relevant genes and/or an exhaustive phenotypic characterisation. This, in turn, suggests that the use of expression data to infer the molecular basis of phenotypes is more challenging than hitherto appreciated. A major advantage of this study is the use of individual-level gene expression data from a large number of subjects, including individuals reared in a shared social environment but exhibiting very different phenotypic responses to perturbation. By sequencing individuals rather than pools, we were able to match each gene expression profile to high-resolution phenotypic data that captures the scale of naturally-occurring variation in features such as age, ovarian development and dominance behavior. This resolution allows us to address questions that are otherwise inaccessible in gene expression analyses. For example, we have been able to show that caste identity, but not the residuals of caste identity on age, are significantly predictive of individuals’ change in transcriptomic caste identity following queen loss. Our discovery of colony-wide responses to queen loss suggests that this social perturbation provokes a significant reaction even from individuals that have little hope of attaining the vacant reproductive role. This is a surprising finding given that P. dominula is thought to express a ‘conventional’ gerontocratic mechanism of dominance and queen succession that mitigates the need for costly intragroup conflicts over the identity of the replacement queen (Pardi 1948 ; Tsuji & Tsuji 2005 ; Monnin et al 2009 ), which should greatly reduce the need for young, low-ranking workers to respond to queen loss (Taylor et al 2020 ). The gene expression responses of lower-ranked workers to queen removal might plausibly represent a form of safeguard against queen loss: if queen loss sometimes occurs multiple times in quick succession or is frequently associated with a general decimation of the nest population (i.e. through predation), there might be kin-selected benefits of a colony-wide ‘de-differentiation’ of individuals that facilitates a quicker succession process. Replacement queens in our colonies did not have access to males and therefore remained unmated even after queen succession. This may partially explain the fact that even individuals with fully developed ovaries and very high dominance ratings did not transition to a fully queen-like gene expression profile, as mating can induce significant gene expression changes in insects (e.g. Gomulski et al 2012 ; Zhou et al 2014 ). The lack of immediate mating opportunities for new queens in our experiment is not necessarily unrealistic, however: unmated Hymenopteran females can lay unfertilized eggs, which develop as males. Moreover, in naturally-occurring early P. dominula nests, replacement queens may be established a month or more before they are mated (Strassmann et al 2004 ), presumably due to a scarcity of early males. The unmated replacement queens analysed here are therefore representative of those that would be present on wild nests shortly after queen loss. To our knowledge, this study represents the first application of a support vector classification approach to behavior-associated transcriptomic data. Using this approach we identified a large group of genes as differing meaningfully between Polistes castes—over 10% of annotated genes, a much larger set than those found using standard analytic approaches in this study and others (e.g. Ferreira et al 2013 ; Toth et al 2014 ; Patalano et al 2015 ; Geffre et al 2017 ). Standard approaches using packages such as edgeR (Robinson et al 2010 ), DESeq2 (Love et al 2014 ), or NOISeq (Tarazona et al 2011 ), which typically include information-sharing between genes and relatively strict fold change cutoffs, allow differential expression to be assessed with a high degree of confidence at the level of individual genes. However, because such approaches assess genes individually and do not reduce the dimensionality of the samples’ transcriptional profiles, they are not well-suited to the tracking of subtle but consistent broad-scale changes. This is especially true for gene expression data that are noisy and heterogenous, or those that are dominated by many genes of small effect (Huang et al 2018 ). The molecular bases of caste in simple social insect societies represent such a class of data. Using a machine learning approach, here we have undertaken a detailed analysis of the relationship between gene expression and socially-mediated phenotypic plasticity, revealing broad-scale changes in caste-associated gene expression profiles following a major social disruption. Our results reveal a hitherto unrecognized capacity for large scale disruption to caste-biased gene expression profiles even in the absence of apparent changes in caste phenotype, a disconnect that undermines simplistic models of the relationship between transcriptome and phenotype. Future studies should continue to marry detailed phenotypic and gene expression data in order to assess the prevalence and provenance of such discontinuities. Supplementary Methods Gene expression quantification Brain tissue was extracted from the heads of individual samples and RNA was extracted using the RNeasy Mini Kit (Qiagen) according to manufacturer’s instructions. Library preparation was performed by Novogene Co. followed by sequencing on an Illumina HiSeq 2000 platform with 150-base pair paired-end reads. Reads were filtered with SortMeRNA (Kopylova et al 2012) using default options to remove ribosomal sequences. Trimmomatic (Bolger et al 2014) was then used to perform quality trimming. First, we trimmed adapter sequences and leading and trailing bases with low phred scores (<3). We then used the MAXINFO option with target length 36 and strictness 0.7 to trim low-quality sequences from the remaining reads. Reads were next mapped to 11313 transcripts from the P. dominula genome annotation 1.0 (Standage et al. 2016) using STAR (Dobin et al 2013) with default options. All 101 samples produced >85% uniquely mapped reads. Reads were assembled into transcripts using StringTie2 (Kovaka et al 2019) before being passed on as raw counts to downstream analyses. Prior to downstream analysis, transcripts were filtered to remove any gene which did not have >20 counts across all samples in at least one of the experimental groups (queens, control workers, and day three and day twelve post-manipulation workers). Following this filtering, 10734/11313 (94.9%) of transcripts remained. Support vector machine classification Support vector classification was performed in R using the e1071 package (Meyer et al 2017) using the gene expression profiles for all available queens (coded with a value of 1; n=26) and control workers (coded with a value of 0; n=12). Classifiers were assessed via their 3-fold cross-validation error rates; classifiers with lower classification error were considered to be superior. Initially, we tested classifiers using a variety of kernel functions (radial, linear, sigmoid and polynomial) combined with grid searches across a wide range of kernel parameters. A radial kernel with γ=10 -6 and cost parameter C=2 5 was found to produce the lowest error rate of all combinations, so all subsequent classifiers were fit using this kernel and a more focused grid search in a parameter space of 2 4 <C<2 6 and 10 -5 <γ<10 -7 in order to minimise the processing power necessary to perform feature selection. For feature selection, we took the classifier fitted with all genes, and iteratively performed the following process: (1) The 3-fold cross-validation error of the model was calculated twenty times using randomly-assigned bins, and the mean of the resulting errors was recorded as the true validation error of the classifier; (2) The feature weights of all genes were calculated by taking the matrix product of that classifier’s coefficients with its support vectors; (3) The gene with the smallest absolute weight in the model was dropped; (4) A new classifier was calculated using the remaining set of genes. This process was repeated until just 100 genes remained, and the optimal support vector classifier was then taken as that for which the cross-validation error reached its minimum. Differential expression analysis Differential expression analyses were performed in R using the DESeq2 package (Love et al 2014). DESeq2 was run on all groups and contrasts were then calculated for each pair of groups. Differential expression was calculated relative to a baseline fold change of 1.5, i.e. p -values refer to the probability that absolute change between two groups was greater than 50%. Genes were considered differentially expressed between conditions if p <0.05 after false discovery rate correction according to the Benjamini-Hochberg procedure. Gene co-expression network analysis Weighted gene co-expression network analysis was performed in R using the WGCNA package (Langfelder & Horvath 2008). As WGCNA is particularly sensitive to genes with low expression, data were first subjected to a second round of filtering in which genes that had 90% of samples were removed. Consensus gene modules across all samples were then constructed using a soft-threshold power of 9. Initially, 26 gene modules were identified. Modules whose eigengene correlation was >0.75 were subsequently merged, after which 22 consensus modules remained. Finally, the Pearson correlation of each module with each phenotypic trait within each group (queens, control workers and individuals from queen removal nests) was calculated and subjected to Benjamini-Hochberg FDR correction. Gene ontology (GO) enrichment analysis In order to perform GO enrichment analysis, we first used OrthoFinder (Emms & Kelly 2019) to identify orthologues for each P. dominula gene in D. melanogaster , a model species for which GO annotations are much more complete. GO annotations for each D. melanogaster gene were acquired from BioMart (Smedley et al 2015) and each P. dominula gene was then assigned GO terms permissively, i.e. a given P. dominula gene was assigned a GO term if that term appeared as an annotation to any of its orthologues. GO enrichment analysis was then performed in R via the topGO package (Alexa & Rahnenfuhrer 2009) using TopGO’s weight01 algorithm and Fisher’s exact test to identify GO terms that were significantly overrepresented ( p <0.01) in a focal set of genes against a background consisting of all genes that appeared in the relevant analysis. Declarations Author contributions BAT, SS and MR conceived of and designed the experiment. BAT and AC carried out the initial experiment. BAT performed behavioural analyses, dissections and RNA extractions, and wrote the initial draft of the manuscript. BAT conducted the analyses with advice from CDRW, MR and SS. BAT and CDRW designed the figures. All authors read and contributed to the manuscript. Acknowledgements We would like to thank F. Cappa, I. Pepiciello, L. Dapporto and others at the Dipartimento di Biologia, Università di Firenze for assistance with field work, and E. Favreau and M. Bentley for their assistance with analyses. Special thanks to R. Cervo, without whose generous support this study would not have been possible. This work was supported by the Natural Environment Research Council (grant code NE/L002485/1). Data availability Sequencing data associated with this paper have been deposited in the NCBI Gene Expression Omnibus (www.ncbi.nlm.nih.gov/geo) under accession number GSE153532 with reviewer access token adctccsahjgjbqt . References Abeel, T., Helleputte, T., Van de Peer, Y., Dupont, P., & Saeys, Y. Robust biomarker identification for cancer diagnosis with ensemble feature selection methods. Bioinformatics 26 , 392–398 (2010). Alexa, A., & Rahnenführer, J. Gene set enrichment analysis with topGO. Bioconductor Improv 27 (2009). Bolger, A. M., Lohse, M., & Usadel, B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30 , 2114–2120 (2014). Casas, L., Saborido-Rey, F., Ryu, T., Michell, C., Ravasi, T., & Irigoien, X. Sex change in clownfish: molecular insights from transcriptome analysis. Sci. Rep. 6 , 35461. (2016). Chakravarty, P., Cozzi, G., Ozgul, A., & Aminian, K. A novel biomechanical approach for animal behaviour recognition using accelerometers. Methods Ecol. Evol. 10 , 802–814 (2019). Chandra, V., Fetter-Pruneda, I., Oxley, P. R., Ritger, A. L., McKenzie, S. K., Libbrecht, R., & Kronauer, D. J. Social regulation of insulin signaling and the evolution of eusociality in ants. Science 361 , 398–402 (2018). Chapman, B. B., Brönmark, C., Nilsson, J. Å., & Hansson, L. A. The ecology and evolution of partial migration. Oikos 120 , 1764–1775 (2011). Cini, A., Meconcelli, S., & Cervo, R. Ovarian indexes as indicators of reproductive investment and egg-laying activity in social insects: a comparison among methods. Insectes Soc. 60 , 393–402 (2013). Cullen, D. A., Cease, A. J., Latchininsky, A. V., Ayali, A., Berry, K., Buhl, J., … Rogers, S. M. From molecules to management: mechanisms and consequences of locust phase polyphenism. Adv. Insect Physiol 53 , 167–285 (2017). Dapporto, L., Sledge, F. M., & Turillazzi, S. Dynamics of cuticular chemical profiles of Polistes dominulus workers in orphaned nests (Hymenoptera, Vespidae). J. Insect Physiol. 51 , 969–973 (2005). Dobin, A., Davis, C. A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., … Gingeras, T. R. STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29 , 15–21 (2013). Elo, A. E. (1978). The rating of chess players, past and present . New York: Arco Pub. Emms, D. M., & Kelly, S. OrthoFinder: phylogenetic orthology inference for comparative genomics. Genome Biol. 20 , 1–14 (2019). Ferreira, P. G., Patalano, S., Chauhan, R., Ffrench-Constant, R., Gabaldon, T., Guigo, R., & Sumner, S. Transcriptome analyses of primitively eusocial wasps reveal novel insights into the evolution of sociality and the origin of alternative phenotypes. Genome Biol. 14 , R20 (2013). Fowler-Finn, K. D., & Rodríguez, R. L. The evolution of experience‐mediated plasticity in mate preferences. J. Evol. Biol. 25 , 1855–1863 (2012). Geffre, A. C., Liu, R., Manfredini, F., Beani, L., Kathirithamby, J., Grozinger, C. M., & Toth, A. L. Transcriptomics of an extended phenotype: parasite manipulation of wasp social behaviour shifts expression of caste-related genes. Proc. R. Soc. B 284 , 20170029 (2017). Gomulski, L. M., Dimopoulos, G., Xi, Z., Scolari, F., Gabrieli, P., Siciliano, P., … Gasperi, G. Transcriptome profiling of sexual maturation and mating in the Mediterranean fruit fly, Ceratitis capitata . PLoS One 7 , (2012). Hadfield, J. D., Nutall, A., Osorio, D., & Owens, I. P. F. Testing the phenotypic gambit: phenotypic, genetic and environmental correlations of colour. J. Evol. Biol. 20 , 549–557 (2007). Heyes, C. Blackboxing: social learning strategies and cultural evolution. Philos. Trans. R. Soc. B 371 , 20150369 (2016). Horiguchi, R., Nozu, R., Hirai, T., Kobayashi, Y., Nagahama, Y., & Nakamura, M. Characterization of gonadal soma-derived factor expression during sex change in the protogynous wrasse, Halichoeres trimaculatus . Dev. Dynam. 242 , 388–399 (2013). Huang, S., Cai, N., Pacheco, P. P., Narrandes, S., Wang, Y., & Xu, W. Applications of support vector machine (SVM) learning in cancer genomics. Cancer Genom. Proteom. 15 , 41–51 (2018). Jandt, J. M., Tibbetts, E. A., & Toth, A. L. Polistes paper wasps: a model genus for the study of social dominance hierarchies. Insectes Soc. 61 , 11–27 (2014). Kopylova, E., Noé, L., & Touzet, H. SortMeRNA: fast and accurate filtering of ribosomal RNAs in metatranscriptomic data. Bioinformatics 28 , 3211–3217 (2012). Kovaka, S., Zimin, A. V., Pertea, G. M., Razaghi, R., Salzberg, S. L., & Pertea, M. Transcriptome assembly from long-read RNA-seq alignments with StringTie2. Genome Biol. 20 , 1–13 (2019). Langfelder, P., & Horvath, S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics 9 , 559 (2008). Libbrecht, R., Oxley, P. R., & Kronauer, D. J. Clonal raider ant brain transcriptomics identifies candidate molecular mechanisms for reproductive division of labor. BMC Biol. 16 , 89 (2018). Lo, N., Simpson, S. J., & Sword, G. A. Epigenetics and developmental plasticity in orthopteroid insects. Curr. Opin. Insect. Sci. 25 , 25–34 (2018). Love, M., Anders, S., & Huber, W. Differential analysis of count data–the DESeq2 package. Genome Biol. 15 , 10–1186 (2014). Manfredini, F., Brown, M. J., & Toth, A. L. Candidate genes for cooperation and aggression in the social wasp Polistes dominula . J. Comp. Physiol. A 204 , 449–463 (2018). Matson, C. K., Murphy, M. W., Sarver, A. L., Griswold, M. D., Bardwell, V. J., & Zarkower, D. DMRT1 prevents female reprogramming in the postnatal mammalian testis. Nature 476 , 101–104 (2011). Meyer, D., Dimitriadou, E., Hornik, K., Weingessel, A. & Leisch, F. e1071: Misc Functions of the Department of Statistics, Probability Theory Group (Formerly: E1071), TU Wien. R package version 1.6-8. https://CRAN.R-project.org/package=e1071 . (2017). Miyano, S. Colony development, worker behavior and male production in orphan colonies of a Japanese paper wasp, Polistes chinensis antennalis Pérez. Popul. Ecol. 68 , 347–361 (1986). Miyano, S. Worker reproduction and related behavior in orphan colonies of a Japanese paper wasp, Polistes jadwigae (Hymenoptera, Vespidae). J. Ethol. 9 , 135–146 (1991). Monnin, T., Cini, A., Lecat, V., Fédérici, P., & Doums, C. No actual conflict over colony inheritance despite high potential conflict in the social wasp Polistes dominulus . Proc. R. Soc. B 276 , 1593–1601 (2009). Neumann, C., Duboscq, J., Dubuc, C., Ginting, A., Irwan, A. M., Agil, M., … Engelhardt, A. Assessing dominance hierarchies: validation and advantages of progressive evaluation with Elo-rating. Anim. Behav. 82 , 911–921 (2011). Pardi, L. Dominance order in Polistes wasps. Physiol. Zool. 21 , 1–13 (1948). Patalano, S., Vlasova, A., Wyatt, C., Ewels, P., Camara, F., Ferreira, P. G., … Sumner, S. Molecular signatures of plastic phenotypes in two eusocial insect species with simple societies. Proc. Natl. Acad. Sci. U.S.A. 112 , 13970–13975 (2015). Rast, W., Kimmig, S. E., Giese, L., & Berger, A. Machine learning goes wild: Using data from captive individuals to infer wildlife behaviours. PLoS One 15 , e0227317 (2020). Réale, D., Dingemanse, N. J., Kazem, A. J., & Wright, J. Evolutionary and ecological approaches to the study of personality. Philos. Trans. R. Soc. B. 365 , 3947–3946 (2010). Rehan, S. M., Glastad, K. M., Steffen, M. A., Fay, C. R., Hunt, B. G., & Toth, A. L. Conserved genes underlie phenotypic plasticity in an incipiently social bee. Genome Biol. Evol. 10 , 2749–2758 (2018). Rittschof, C. C., & Robinson, G. E. Genomics: moving behavioural ecology beyond the phenotypic gambit. Anim. Behav. 92 , 263–270 (2014). Robinson, M. D., McCarthy, D. J., & Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics 26 , 139–140 (2010). Rubin, H. The phenotypic gambit: selective pressures and ESS methodology in evolutionary game theory. Biol. Philos. 31 , 551–569 (2016). Schwander, T., Lo, N., Beekman, M., Oldroyd, B. P., & Keller, L. Nature versus nurture in social insect caste differentiation. Trends Ecol. Evol. 25 , 275–282 (2010). Segal, N. H., Pavlidis, P., Antonescu, C. R., Maki, R. G., Noble, W. S., DeSantis, D., … Cordon-Cardo, C. Classification and subtype prediction of adult soft tissue sarcoma by functional genomics. Am. J. Pathol. 163 , 691–700 (2003). Shell, W. A., & Rehan, S. M. Behavioral and genetic mechanisms of social evolution: insights from incipiently and facultatively social bees. Apidologie 49 , 13–30 (2018). Simola, D. F., Graham, R. J., Brady, C. M., Enzmann, B. L., Desplan, C., Ray, A., … Berger, S. L. Epigenetic (re) programming of caste-specific behavior in the ant Camponotus floridanus . Science 351 , aac6633 (2016). Smedley, D., Haider, S., Durinck, S., Pandini, L., Provero, P., Allen, J., … Bardou, P. The BioMart community portal: an innovative alternative to large, centralized data repositories. Nucleic Acids Res. 43 , W589-W598 (2015). Standage, D. S., Berens, A. J., Glastad, K. M., Severin, A. J., Brendel, V. P., & Toth, A. L. Genome, transcriptome and methylome sequencing of a primitively eusocial wasp reveal a greatly reduced DNA methylation system in a social insect. Mol. Ecol. 25 , 1769–1784 (2016). Starks, P. T., Turillazzi, S., & West-Eberhard, M. J. Polistes paper wasps: emergence of a model genus. In Ann. Zool. Fenn. 43 , 385–386 (2006). Strassmann, J. E., Fortunato, A., Cervo, R., Turillazzi, S., Damon, J. M., & Queller, D. C. The cost of queen loss in the social wasp Polistes dominulus (Hymenoptera: Vespidae). J. Kansas Entomol. Soc. 77 , 343–355 (2004). Tarazona, S., García, F., Ferrer, A., Dopazo, J., & Conesa, A. NOIseq: a RNA-seq differential expression method robust for sequencing depth biases. EMBnet.journal 17 , 18–19. (2011). Taylor, B. A., Cini, A., Cervo, R., Reuter, M., & Sumner, S. Queen succession conflict in the paper wasp Polistes dominula is mitigated by age-based convention. Behav. Ecol. , araa045 (2020). Tibbetts, E. A., & Huang, Z. Y. The challenge hypothesis in an insect: juvenile hormone increases during reproductive conflict following queen loss in Polistes wasps. Am. Nat. 176 , 123–130 (2010). Todd, E. V., Ortega-Recalde, O., Liu, H., Lamm, M. S., Rutherford, K. M., Cross, H., … Godwin, J. R. Stress, novel sex genes, and epigenetic reprogramming orchestrate socially controlled sex change. Sci. Adv. 5 , eaaw7006 (2019). Toth, A. L., Tooker, J. F., Radhakrishnan, S., Minard, R., Henshaw, M. T., & Grozinger, C. M. Shared genes related to aggression, rather than chemical communication, are associated with reproductive dominance in paper wasps ( Polistes metricus ). BMC Genom. 15 , 75 (2014). Tsuji, K., & Tsuji, N. Why is dominance hierarchy age-related in social insects? The relative longevity hypothesis. Behav. Ecol. Sociobiol. 58 , 517–526 (2005). West-Eberhard, M. J. The Social Biology of Polistine Wasps. Misc. Pub. Mus. Zool. Univ. Mich. 140 , 1–101 (1969). Zhou, S., Mackay, T. F., & Anholt, R. R. Transcriptional and epigenetic responses to mating and aging in Drosophila melanogaster . BMC Genom. 15 , 927 (2014). Additional Declarations There is NO Competing Interest. Supplementary Files supplementarytableS1.pdf Supplementary Table S1 supplementarytableS2.pdf Supplementary Table S2 supplementarytableS3.pdf Supplementary Table S3 supplementarytableS4.pdf Supplementary Table S4 supplementaryfigures.pdf Supplementary Figure S1 Cite Share Download PDF Status: Published Journal Publication published 02 Feb, 2021 Read the published version in Nature Communications → Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-43815","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":898873,"identity":"c8346e2c-ebc5-4b64-9a8c-54cc942ae63c","order_by":0,"name":"Benjamin Taylor","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAABLElEQVRIie2RMUvDQBSAn1Quy6nrk2rzFy4cpA0o/gF/xAXXIIUsmaQQaJeErv0ZipA5INQlOmfSunSy0KnUQfCSM9sls+B98Lh3x3137+4BGAx/EKYC1Wwlg1hxDjCu571uhchMyDimS6HWOxRFowzwmnUqQyvNx1/B6BbseI0iegOCdId7BvbJhHKmUbzkVThphiEQ4qIo5NhPs9OEgbPIKRe6wsqAsaMM/QkBF/2pAHL2kvVlYQf3QHmuU94/mfNdK9ZOKRisK+WqVSkp4+oW6jYKqRS/UnSFeUnA+LlUpiQIR6IQlNAl9xKGN4snEuqeP7QK19lkd/689/xQbiMxsGfxR7mPLi7ns/gRdb8sG1IfJdtyWG2gv6vY0pWa3qpJtq17DAaD4V/zA4dRVCkp/jJNAAAAAElFTkSuQmCC","orcid":"","institution":"University College London","correspondingAuthor":true,"submittingAuthor":false,"prefix":"","firstName":"Benjamin","middleName":"","lastName":"Taylor","suffix":""},{"id":898874,"identity":"199966f1-ae54-4544-a345-942b6dc3f4d7","order_by":1,"name":"Alessandro Cini","email":"","orcid":"https://orcid.org/0000-0003-0355-2188","institution":"Centre for Biodiversity and Environment Research, University College London","correspondingAuthor":false,"submittingAuthor":false,"prefix":"","firstName":"Alessandro","middleName":"","lastName":"Cini","suffix":""},{"id":898875,"identity":"f1c38c21-dfa8-4704-a021-905eae2974b2","order_by":2,"name":"Christopher Wyatt","email":"","orcid":"","institution":"Centre for Genomic Regulation (CRG), The Barcelona Institute of Science and Technology, Barcelona","correspondingAuthor":false,"submittingAuthor":false,"prefix":"","firstName":"Christopher","middleName":"","lastName":"Wyatt","suffix":""},{"id":898876,"identity":"154c4c6b-3106-41d8-8ed0-c917e8a04b1e","order_by":3,"name":"Max Reuter","email":"","orcid":"https://orcid.org/0000-0001-9554-0795","institution":"[email protected]","correspondingAuthor":false,"submittingAuthor":false,"prefix":"","firstName":"Max","middleName":"","lastName":"Reuter","suffix":""},{"id":898877,"identity":"408478af-3cdb-447e-ae2a-56355bbe4ae8","order_by":4,"name":"Seirian Sumner","email":"","orcid":"https://orcid.org/0000-0003-0213-2018","institution":"University College London","correspondingAuthor":false,"submittingAuthor":false,"prefix":"","firstName":"Seirian","middleName":"","lastName":"Sumner","suffix":""}],"badges":[],"createdAt":"2020-07-15 15:31:39","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-43815/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-43815/v1","draftVersion":[],"editorialEvents":[{"content":"https://doi.org/10.1038/s41467-021-21095-6","type":"published","date":"2021-02-03T02:42:17+00:00"}],"editorialNote":"","failedWorkflow":false,"files":[{"id":1696789,"identity":"fd2623b1-1dc8-4e73-adac-762520fdf016","added_by":"auto","created_at":"2020-07-27 18:44:00","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":33368,"visible":true,"origin":"","legend":"Classifier estimates generated by the lowest-error SVM using 1992 caste-informative genes for (A) queens and workers from control colonies and (B) individuals from experimental colonies following queen removal. Mean and standard deviation for each group are shown in black.","description":"","filename":"Onlinefloatimage2.Png","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/Onlinefloatimage2.Png"},{"id":1696790,"identity":"4fea516e-9b48-4e7f-b67c-39007fcc6fcc","added_by":"auto","created_at":"2020-07-27 18:44:00","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":102730,"visible":true,"origin":"","legend":"(A-C) Scatterplots of SVM classifier estimates for individuals from post-removal colonies plotted against (A) phenotypic queenness, (B) age, and (C) the residuals of queenness on age respectively. (D-F) Module-trait correlations for consensus modules identified by WGCNA in queen removal individuals, control workers and queens respectively. Within each cell, the Pearson correlation coefficient is given above and the FDR-adjusted p-value is given below in parentheses. (G-H) Enriched gene ontology terms at p\u003c0.01 for gene module 11 and gene module 8.","description":"","filename":"Onlinefloatimage3.Png","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/Onlinefloatimage3.Png"},{"id":1696791,"identity":"582c1fd9-1ed1-4410-9759-c27d73bc4fb6","added_by":"auto","created_at":"2020-07-27 18:44:00","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":78606,"visible":true,"origin":"","legend":"Box 1. Summarized experimental design describing the generation of phenotypic data pertaining to individual-level responses to queen removal in Polistes dominula (Taylor et al 2020) used in this study. (A) Early-season nests were transferred to the lab prior to the eclosion of workers and subordinate foundresses were removed. After at least 4 workers had eclosed, nests were filmed for three days and then assigned to either control or queen removal treatments. Following treatment (sham or queen removal) nests were filmed for a further three or twelve days and then all individuals were dissected to generate ovarian development indices following Cini et al (Cini et al. 2013). Footage of nests was used to generate dominance indices in the form of Elo ratings (Elo 1978; Neumann et al. 2011). (B) Ovarian and dominance indices from queens and control workers were used to produce a logistic regression model for caste classification (0=worker, 1=queen). Data from individuals on queenless post-removal nests were then passed through this model to fit caste estimates and thereby identify individuals with high ‘queenness’, i.e. those that exhibited strongly queen-like phenotypes following queen removal and thus represented possible replacement queens. Of the individuals for which data are shown here, 27 queens, 12 workers from control nests, and 62 individuals from queen removal nests were subsequently sequenced to generate the data discussed in the present study.","description":"","filename":"Onlinefloatimage1.Png","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/Onlinefloatimage1.Png"},{"id":13559843,"identity":"ff6bcf68-b453-4628-ae56-a33635c67dfe","added_by":"auto","created_at":"2021-09-17 03:02:18","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":1130743,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/43c6ebe8-a429-469d-8062-cadd125ff049.pdf"},{"id":1696793,"identity":"6ecfbcf5-efdd-4782-b364-0e9bbf516270","added_by":"auto","created_at":"2020-07-27 18:44:03","extension":"pdf","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":6373197,"visible":true,"origin":"","legend":"Supplementary Table S1","description":"","filename":"supplementarytableS1.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/supplementarytableS1.pdf"},{"id":1696794,"identity":"3e4450a6-b1c6-4a64-b854-881f69cf8653","added_by":"auto","created_at":"2020-07-27 18:44:03","extension":"pdf","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":791619,"visible":true,"origin":"","legend":"Supplementary Table S2","description":"","filename":"supplementarytableS2.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/supplementarytableS2.pdf"},{"id":1696795,"identity":"ad582d0d-3559-40b9-a1f5-ba285c195836","added_by":"auto","created_at":"2020-07-27 18:44:03","extension":"pdf","order_by":3,"title":"","display":"","copyAsset":false,"role":"supplement","size":223410,"visible":true,"origin":"","legend":"Supplementary Table S3","description":"","filename":"supplementarytableS3.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/supplementarytableS3.pdf"},{"id":1696796,"identity":"c39325b5-2808-4219-8d02-c7b2823d3af8","added_by":"auto","created_at":"2020-07-27 18:44:04","extension":"pdf","order_by":4,"title":"","display":"","copyAsset":false,"role":"supplement","size":158551,"visible":true,"origin":"","legend":"Supplementary Table S4","description":"","filename":"supplementarytableS4.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/supplementarytableS4.pdf"},{"id":1696797,"identity":"5668df78-e303-47f2-9c12-e31c2a045dec","added_by":"auto","created_at":"2020-07-27 18:44:04","extension":"pdf","order_by":5,"title":"","display":"","copyAsset":false,"role":"supplement","size":310963,"visible":true,"origin":"","legend":"Supplementary Figure S1","description":"","filename":"supplementaryfigures.pdf","url":"https://assets-eu.researchsquare.com/files/rs-43815/v1/supplementaryfigures.pdf"}],"financialInterests":"There is \u003cb\u003eNO\u003c/b\u003e Competing Interest.","formattedTitle":"The molecular basis of socially-mediated phenotypic plasticity in a eusocial paper wasp","fulltext":[{"header":"Introduction","content":" \u003cp\u003eThe relationship between gene expression and external phenotype is complex and unresolved. Much research in behavioural and evolutionary ecology is based on the implicit assumption that phenotypic traits can be modelled as though they directly reflect gene expression patterns, and that evolutionary trajectories can therefore be studied while remaining agnostic with regard to the underlying molecular mechanisms (Hadfield et al \u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e2007\u003c/span\u003e; Rittschof \u0026amp; Robinson \u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e2014\u003c/span\u003e; Rubin \u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e2016\u003c/span\u003e). This \u0026lsquo;phenotypic gambit\u0026rsquo; has proven a useful rule of thumb, permitting the establishment of a rich body of literature surrounding the evolution of complex traits, despite a lack of data relating to the genetic basis of these traits (e.g. R\u0026eacute;ale et al \u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e2010\u003c/span\u003e; Chapman et al \u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e2011\u003c/span\u003e; Fowler-Finn \u0026amp; Rodr\u0026iacute;guez \u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e2012\u003c/span\u003e). In the past decade, however, advances in the affordability of \u0026lsquo;omic\u0026rsquo; data and availability of powerful bioinformatic methods have greatly enhanced our ability to assess the assumptions made by the phenotypic gambit (Rittschof \u0026amp; Robinson \u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e2014\u003c/span\u003e; Heyes \u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e2016\u003c/span\u003e). The time is right to disentangle the molecular bases of complex phenotypic traits.\u003c/p\u003e \u003cp\u003ePhenotypic plasticity, the ability of an individual to effect phenotypic changes in response to external cues, is an ideal phenomenon with which to study the relationship between gene expression and phenotype because it involves the production of multiple phenotypes without gene sequence changes. Of particular value are species in which adult individuals can be experimentally induced to transition between distinct, measurable phenotypes. By comparing the gene expression profiles of groups of individuals that differ in their phenotypes as a result of plasticity, it is possible to isolate phenotypic effects of gene expression while controlling for genetic differences. Using this approach, significant progress has been made in unravelling the molecular underpinnings of sequential sex changes in hermaphroditic fish (Horiguchi et al \u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e2013\u003c/span\u003e; Casas et al \u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e2016\u003c/span\u003e; Todd et al \u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e2019\u003c/span\u003e), the distinct gregarious social phenotype of desert locusts (Cullen et al \u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e2017\u003c/span\u003e; Lo et al \u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e2018\u003c/span\u003e), and the reproductive castes of social insects (Simola et al \u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e2016\u003c/span\u003e; Libbrecht et al \u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e2018\u003c/span\u003e; Rehan et al \u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e2018\u003c/span\u003e; Shell \u0026amp; Rehan \u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e2018\u003c/span\u003e). Such studies typically rely on comparisons between groups of individuals with well-differentiated phenotypes, however. As a result, little is known about more subtle effects during the transition from one morph to the other, and the relationship between expression patterns and phenotypic traits at the individual level.\u003c/p\u003e \u003cp\u003eThe reproductive castes found in the colonies of social insects provide excellent model systems for determining the extent to which fine-scale changes in phenotype are reflected at the molecular level. With a few exceptions (reviewed in Schwander et al \u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e2010\u003c/span\u003e), the distinct queen and worker phenotypes found in such colonies are plastically determined either during development or in adulthood. Workers in some species can be experimentally induced to transition to a reproductive role in response to the removal of a colony\u0026rsquo;s queen (e.g. Strassmann et al \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e2004\u003c/span\u003e; Tibbetts \u0026amp; Huang \u003cspan citationid=\"CR54\" class=\"CitationRef\"\u003e2010\u003c/span\u003e) or as a result of exposure to varying levels of brood (e.g. Chandra et al \u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e2018\u003c/span\u003e; Libbrecht et al \u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e2018\u003c/span\u003e), allowing changes in the behavioural, physiological and molecular traits that define caste identity to be tracked. An additional benefit\u0026mdash;and challenge\u0026mdash;of studying social insect colonies is that they involve complex social structures. Such interactions can be hard to study, but offer the opportunity to assess the effects of social interactions upon phenotypes and transcriptomes.\u003c/p\u003e \u003cp\u003eIn the present study, we explore the relationship between transcriptome and phenotype in the context of plastic caste expression in the European paper wasp \u003cem\u003ePolistes dominula\u003c/em\u003e (Christ 1791), which is a model organism often used in studies of social insect behaviour (reviewed in Starks \u0026amp; Turillazzi 2006; Jandt et al \u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e2014\u003c/span\u003e) and, more recently, for analyses of caste gene expression (e.g. Patalano et al \u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e2015\u003c/span\u003e; Standage et al \u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e2016\u003c/span\u003e; Geffre et al \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e2017\u003c/span\u003e; Manfredini et al \u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e2018\u003c/span\u003e). In this species, removing the established queen from a single-foundress colony induces a queen succession process in which one (or very few) workers transition to a queen phenotype, with age playing a key role in predicting which individual will do so (Pardi \u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e1948\u003c/span\u003e; Tsuji \u0026amp;Tsuji \u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e2005\u003c/span\u003e): almost invariably, the new queen is one of the oldest individuals, and there is little conflict over succession (Strassmann et al \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e2004\u003c/span\u003e; Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e). In a recent paper (Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e) we followed responses to queen removal in \u003cem\u003eP. dominula\u003c/em\u003e on a fine scale by measuring individual-level behavioural and physiological traits and generating a univariate measure of caste identity (\u0026lsquo;queenness\u0026rsquo;) to describe individuals\u0026rsquo; phenotypic profiles (\u003cb\u003eBox 1\u003c/b\u003e). We found that phenotypic responses were limited to the subset of individuals that transitioned to the queen role, while other individuals exhibited little or no measurable change in measured traits. This system, in which individuals within a controlled environment vary strongly and predictably in their phenotypic response to a shared stimulus, affords an excellent opportunity to track the relationship between gene expression and phenotypic expression at the individual level.\u003c/p\u003e \u003cp\u003eTo do so, here we complement the individual-level phenotypic data of Taylor et al. (\u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e) by assaying the brain gene expression profiles of 107 individuals for which fine-scale behavioural and ovarian data have been obtained, including queens and workers from stable colonies and individuals from colonies that had their queens experimentally removed. We then assess the relationship between the transcriptome and the phenotypic changes induced by queen removal. To match the univariate measure of caste identity that we applied to our phenotypic data (\u003cb\u003eBox 1\u003c/b\u003e), we use a support vector classification approach to generate a univariate measure of caste-associated gene expression. Support vector classification is a powerful tool with which to transform complex patterns in multidimensional data into a continuous classification score, allowing the detection of subtle, widespread signals of differential expression between phenotypic states that are likely to be missed in conventional differential expression analyses. Support vector machines (SVMs) have become a key tool in the early identification of phenotypically indistinguishable cancer subtypes (Segal et al \u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e2003\u003c/span\u003e; Abeel et al \u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e2010\u003c/span\u003e; Huang et al \u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e2018\u003c/span\u003e) and their value has recently been demonstrated in animal behaviour studies (e.g. Chakravarty et al \u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e2019\u003c/span\u003e; Rast et al \u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e2020\u003c/span\u003e). Given their sensitivity, SVMs should be ideally suited to quantify gene expression variation across the spectrum between differentiated worker and queen roles, but to our knowledge this is the first time SVMs have been used to interrogate the molecular basis of behavioural traits.\u003c/p\u003e \u003cp\u003eApplying an SVM approach along with standard differential expression and gene co-expression analyses, we show that brain gene expression responses to queen removal in \u003cem\u003eP. dominula\u003c/em\u003e include a colony-wide response that does not match that observed at the phenotypic level. Our results indicate that gene expression in \u003cem\u003eP. dominula\u003c/em\u003e colonies reflects both a generalised response to queen loss which is independent of phenotype, and a phenotype-specific response that tracks individuals\u0026rsquo; expression of plastic phenotypic changes. To our knowledge, this study provides the most comprehensive analysis to date of the ways that plastic phenotypes are reflected at the transcriptomic level; our results expose the complexity of the relationship between individual-level gene expression and individual phenotypes and call into question the assumptions of the phenotypic gambit.\u003c/p\u003e"},{"header":"Results","content":" \u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003eSupport vector classification reveals consistent patterns of caste gene expression differentiation involving many genes\u003c/h2\u003e \u003cp\u003eOur support vector approach allowed us to reduce complex variation in brain gene expression data down to a single dimension of caste identity, matching the Bayesian logistic regression model we used in Taylor et al. (\u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e) to condense ovarian and behavioural data into a unidimensional metric of phenotypic caste (\u0026lsquo;queenness\u0026rsquo;; \u003cb\u003eBox 1\u003c/b\u003e). We trained an SVM using gene expression data from 27 queens and 12 workers from stable, unmanipulated colonies. An initial model, based on all 10734 genes annotated in our experiment, achieved a root mean squared validation error of 0.065 in three-fold cross-validation, i.e. a model trained on a random subset of two thirds of workers and queens classified the remaining third within 25.5% of their true values (worker\u0026thinsp;=\u0026thinsp;0, queen\u0026thinsp;=\u0026thinsp;1). By applying a process of iterative feature selection we were able to remove genes that had low weights in the classification and contributed to overfitting (\u003cb\u003eFigure S1A; Table S1\u003c/b\u003e). Doing so allowed us to identify an optimal model containing 1992 genes with a substantially reduced root mean squared classification error of 0.021 (\u003cb\u003eTable S2\u003c/b\u003e). This model classified queens and workers very consistently, with strong separation of queens from workers (Fig.\u0026nbsp;1A). Thirty-four gene ontology (GO) terms were significantly enriched among these 1992 genes, including a number of terms associated with translation such as \u003cem\u003erRNA processing\u003c/em\u003e, \u003cem\u003etRNA aminoacylation\u003c/em\u003e, and \u003cem\u003eribosomal large subunit biogenesis\u003c/em\u003e (\u003cb\u003eFigure S1B\u003c/b\u003e).\u003c/p\u003e \u003cp\u003eThe 1992 genes in our optimised SVM represent a much larger set than that which is found to differentiate queens and workers based on a standard differential expression approach. When we applied a DESeq2 analysis with a 1.5 fold-change threshold to the same set of queens and workers used to train the SVM, we identified just 81 differentially expressed genes (with no associated GO terms), a small number that is typical for similar analyses in \u003cem\u003ePolistes\u003c/em\u003e (e.g. Toth et al \u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e2014\u003c/span\u003e; Patalano et al \u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e2015\u003c/span\u003e; Geffre et al \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e2017\u003c/span\u003e). All but four of these genes were present in the larger SVM set (\u003cb\u003eFigure S1C\u003c/b\u003e; \u003cb\u003eTable S3\u003c/b\u003e), suggesting that the SVM captures almost all the information obtained from these standard methods. Conversely, an SVM trained using just these 81 differentially-expressed genes exhibited a cross-validation error rate of 0.065, no better than the original un-optimised model. Thus, the picture of caste differentiation provided by standard differential expression analysis appears to miss a great number of subtle differences in individuals\u0026rsquo; gene expression profiles that contribute to caste differentiation.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec4\" class=\"Section2\"\u003e \u003ch2\u003eColony-wide brain gene expression responses to queen removal\u003c/h2\u003e \u003cp\u003eFollowing the loss of a queen from a \u003cem\u003ePolistes\u003c/em\u003e colony, typically one or a few individuals undergo a phenotypic transition to become a replacement queen while the rest of the colony members remain workers (Miyano \u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e1991\u003c/span\u003e; Dapporto et al \u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e2005\u003c/span\u003e; Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e). If this response is reflected at the level of gene expression, then we would expect to observe significant changes in the transcriptional profiles of individuals that exhibit a phenotypic transition, but not in those that do not. Our individual-level expression data allowed us to track responses in 62 individuals\u0026rsquo; gene expression profiles following queen removal across 19 nests. We analysed individuals\u0026rsquo; gene expression profiles at three days after queen removal, when queen replacement is ongoing (n\u0026thinsp;=\u0026thinsp;24), and at twelve days after queen removal, when succession is largely settled at the phenotypic level (Strassmann et al \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e2004\u003c/span\u003e; Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e; n\u0026thinsp;=\u0026thinsp;38). Individuals for sequencing were selected to cover a wide range of phenotypes, including those that remained entirely worker-like, those that had transitioned to highly queen-like phenotypes, and those with intermediate phenotypes at the time of sampling (\u003cb\u003eFigure S1D-F\u003c/b\u003e).\u003c/p\u003e \u003cp\u003eAnalysing shifts in SVM estimates of individual wasps allowed us to assess the degree to which these concurred with the changes visible at the phenotypic level. Doing so, we found that the SVM estimates of individuals from post-removal nests were intermediate between those of queens and workers from control nests (Fig.\u0026nbsp;1B), a finding which concurs with the placement of these individuals according to principal component analysis (\u003cb\u003eFigure S1G\u003c/b\u003e). This result is surprising given that the majority of individuals on queen removal nests are phenotypically indistinguishable from workers on control nests in terms of ovarian development and behavioural dominance (Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e). Thus, the large majority of individuals on queen removal nests exhibited perturbation of their caste-associated gene expression, even though only a few of these individuals exhibited responses to queen loss at the level of physiology or behaviour and transitioned from a worker to a queen role.\u003c/p\u003e \u003cp\u003eFurthermore, we found no evidence that the degree of transcriptional perturbation declined over time following queen removal. SVM classification estimates of individuals from queen removal nests did not differ significantly between days three and twelve following queen removal (QR3 mean 0.331\u0026thinsp;\u0026plusmn;\u0026thinsp;0.111; QR12 mean 0.345\u0026thinsp;\u0026plusmn;\u0026thinsp;0.082; Wilcoxon W\u0026thinsp;=\u0026thinsp;369, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.55). The effects of queen removal on individuals\u0026rsquo; gene expression profiles therefore appear to be both widespread and persistent, affecting all individuals in a nest and lasting beyond the point at which a new queen has already become phenotypically established.\u003c/p\u003e \u003cp\u003eInterestingly, the strong colony-wide perturbation following queen loss that this SVM approach identifies would have been entirely missed using a standard differential expression approach: DESeq2 identified just five genes as differentially expressed between control workers and individuals from manipulated nests, with no associated GO enrichment (\u003cb\u003eTable S4\u003c/b\u003e).\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec5\" class=\"Section2\"\u003e \u003ch2\u003eAge and queenness explain variation in individual-level molecular responses to queen loss\u003c/h2\u003e \u003cp\u003eWhile SVM classification indicates that queen removal causes colony-wide perturbation to brain expression profiles, classifier estimates varied substantially between individuals following queen removal, spanning a much greater range of values (0.116\u0026ndash;0.540) than those of queens (0.900\u0026ndash;1.100) or workers (0.057-0.100) from queenright colonies (Fig.\u0026nbsp;1). To better understand this variation, we examined whether the classifier estimates for individuals from manipulated colonies were predicted by those individuals' phenotypic traits\u0026mdash;specifically ovarian development, behavioural dominance, and age (\u003cb\u003eBox 1\u003c/b\u003e; Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e).\u003c/p\u003e \u003cp\u003eOur phenotypic measure of caste identity (queenness) was a significant predictor of expression-based SVM classifier estimates when fitted using a linear model (slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0414\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0114, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;6.0\u0026thinsp;\u0026times;\u0026thinsp;10\u003csup\u003e\u0026minus;\u0026thinsp;4\u003c/sup\u003e; Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eA). The individual components of queenness, ovarian development and dominance were also significant or near-significant predictors of SVM classification individually (ovarian development: slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0426\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0113, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;4.0\u0026thinsp;\u0026times;\u0026thinsp;10\u003csup\u003e\u0026minus;\u0026thinsp;4\u003c/sup\u003e; dominance: slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0231\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0122, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.06). Notably, however, because phenotypic queenness was strongly correlated with age among post-removal individuals (cor\u0026thinsp;=\u0026thinsp;0.4832, p\u0026thinsp;=\u0026thinsp;8.0\u0026thinsp;\u0026times;\u0026thinsp;10\u003csup\u003e\u0026minus;\u0026thinsp;5\u003c/sup\u003e; Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e), age was an equally strong predictor of caste estimates when fitted in a separate linear model (slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0516\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0106, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;9.8\u0026thinsp;\u0026times;\u0026thinsp;10\u003csup\u003e\u0026minus;\u0026thinsp;6\u003c/sup\u003e; Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eB). Thus, the significance of queenness as a predictor of caste estimates might have been an artefact of the fact that both are correlated with age. To test whether queenness had an effect over and above that accounted for by age, we calculated the residuals of phenotypic queenness on age and fitted the caste estimates against these. These residuals were not significantly predictive of SVM classification (slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0202\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0123, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.11; Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eC), nor were the residuals of ovarian development alone on age (slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0223\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0123, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.07) or the residuals of dominance alone (slope\u0026thinsp;\u0026plusmn;\u0026thinsp;SE\u0026thinsp;=\u0026thinsp;0.0114\u0026thinsp;\u0026plusmn;\u0026thinsp;0.0126, \u003cem\u003ep\u003c/em\u003e\u0026thinsp;=\u0026thinsp;0.37). Age is thus the strongest determinant both of individuals\u0026rsquo; caste phenotypes and of their caste-associated gene expression.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eIn order to further discern the factors that shaped individuals\u0026rsquo; gene expression profiles following queen removal, we performed weighted gene co-expression network analysis (WGCNA) using the full set of annotated genes. We generated 22 consensus modules across all samples, and then determined which traits significantly predicted the expression of a given module within each group (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eD-F). In doing so, we identified three modules of particular interest.\u003c/p\u003e \u003cp\u003eModule 11 (\u003cem\u003eRed\u003c/em\u003e) consists of 519 genes associated with 14 enriched GO terms, including a number of terms associated with molecular binding (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eG). This module is notable in that its expression in individuals from queen removal colonies is positively associated with SVM classifications but not with phenotypic correlates of caste. A second module, Module 2 (\u003cem\u003eLightyellow\u003c/em\u003e), consists of 102 genes and appears to be associated with many caste-related traits in the queen removal condition, being negatively correlated with queenness, age, SVM classification and (more weakly) with ovarian development among queen removal individuals. We also found weak evidence that this module is negatively associated with age among workers from control nests. Unexpectedly, despite the seeming importance of Module 2 in predicting caste identity among workers, not a single GO term was enriched among the genes in the module at p\u0026thinsp;\u0026lt;\u0026thinsp;0.01. Finally, Module 8 (\u003cem\u003eYellow\u003c/em\u003e), consists of 614 genes and has a strongly negative correlation with age among both control workers and queen removal individuals. This module is associated with 21 GO terms, including a number of terms associated with molecular binding, such as protein binding, DNA binding, RNA binding and chromatin binding (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e2\u003c/span\u003eH). The latter term is particularly notable: regulation of chromatin accessibility is one of several epigenetic processes that have been implicated in the control of caste expression in social insects (Simola et al \u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e2016\u003c/span\u003e; Wojciechowski et al 2018; Duncan et al 2020). Downregulation of this module with age might plausibly contribute to the increased phenotypic plasticity exhibited by older individuals in \u003cem\u003eP. dominula.\u003c/em\u003e\u003c/p\u003e \u003c/div\u003e "},{"header":"Discussion","content":" \u003cp\u003eThis study provides one of the most comprehensive analyses to date regarding the transcriptomic signatures of plastic phenotypes, and the first to employ an SVM approach to interrogate the molecular basis of behaviour. We have applied this approach to identify caste-specific gene expression profiles in \u003cem\u003eP. dominula\u003c/em\u003e and to explore the relationship between transcriptomic and plastic phenotypic changes following a major social perturbation in paper wasp colonies. We identify a set of nearly 2000 genes that optimally capture gene expression differences between established queens and workers. Using a caste classifier based on these genes, we find that queen removal leads to a colony-wide shift in expression, where the expression profiles of all individuals move towards a state intermediate between those of established queens and workers. Individual variation around this intermediate state is related to age and phenotypic attributes, with older (and therefore more queen-like) individuals showing expression profiles that are closer to that expected for established queens. Our results show that molecular responses to queen removal in \u003cem\u003eP. dominula\u003c/em\u003e consist of both a general colony-wide response independent of phenotypic change and a response that reflects the plastic phenotypic transition.\u003c/p\u003e \u003cp\u003eOur study contributes to the enormous progress in our understanding of the relationship between molecular changes and changes in phenotypic expression that has been made in the past decade, facilitated by the increased availability of 'omic\u0026rsquo; data and complex bioinformatic analyses. Recent studies have started to challenge the view that there is a direct correspondence between transcriptomic states and external phenotypes. Libbrecht et al (\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e2018\u003c/span\u003e), for example, show that gene expression responses associated with a reversible phenotypic change differ qualitatively based on the directionality of the change (from reproductive to non-reproductive or vice versa). Meanwhile, molecular manipulations have revealed a surprising degree of plasticity in canonically implastic traits such as mammalian sex (Matson et al \u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e2011\u003c/span\u003e) or ant castes (Simola et al \u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e2016\u003c/span\u003e). Our results go further, showing a shift in caste-specific brain gene expression profiles among individuals whose phenotypic caste expression remains otherwise apparently unchanged. This shows that the expectation of a close match between expression profiles and phenotypes is excessively simplistic, or at least that detecting such a match requires detailed knowledge of relevant genes and/or an exhaustive phenotypic characterisation. This, in turn, suggests that the use of expression data to infer the molecular basis of phenotypes is more challenging than hitherto appreciated.\u003c/p\u003e \u003cp\u003eA major advantage of this study is the use of individual-level gene expression data from a large number of subjects, including individuals reared in a shared social environment but exhibiting very different phenotypic responses to perturbation. By sequencing individuals rather than pools, we were able to match each gene expression profile to high-resolution phenotypic data that captures the scale of naturally-occurring variation in features such as age, ovarian development and dominance behavior. This resolution allows us to address questions that are otherwise inaccessible in gene expression analyses. For example, we have been able to show that caste identity, but not the residuals of caste identity on age, are significantly predictive of individuals\u0026rsquo; change in transcriptomic caste identity following queen loss.\u003c/p\u003e \u003cp\u003eOur discovery of colony-wide responses to queen loss suggests that this social perturbation provokes a significant reaction even from individuals that have little hope of attaining the vacant reproductive role. This is a surprising finding given that \u003cem\u003eP. dominula\u003c/em\u003e is thought to express a \u0026lsquo;conventional\u0026rsquo; gerontocratic mechanism of dominance and queen succession that mitigates the need for costly intragroup conflicts over the identity of the replacement queen (Pardi \u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e1948\u003c/span\u003e; Tsuji \u0026amp; Tsuji \u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e2005\u003c/span\u003e; Monnin et al \u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e2009\u003c/span\u003e), which should greatly reduce the need for young, low-ranking workers to respond to queen loss (Taylor et al \u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e2020\u003c/span\u003e). The gene expression responses of lower-ranked workers to queen removal might plausibly represent a form of safeguard against queen loss: if queen loss sometimes occurs multiple times in quick succession or is frequently associated with a general decimation of the nest population (i.e. through predation), there might be kin-selected benefits of a colony-wide \u0026lsquo;de-differentiation\u0026rsquo; of individuals that facilitates a quicker succession process.\u003c/p\u003e \u003cp\u003eReplacement queens in our colonies did not have access to males and therefore remained unmated even after queen succession. This may partially explain the fact that even individuals with fully developed ovaries and very high dominance ratings did not transition to a fully queen-like gene expression profile, as mating can induce significant gene expression changes in insects (e.g. Gomulski et al \u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e2012\u003c/span\u003e; Zhou et al \u003cspan citationid=\"CR59\" class=\"CitationRef\"\u003e2014\u003c/span\u003e). The lack of immediate mating opportunities for new queens in our experiment is not necessarily unrealistic, however: unmated Hymenopteran females can lay unfertilized eggs, which develop as males. Moreover, in naturally-occurring early \u003cem\u003eP. dominula\u003c/em\u003e nests, replacement queens may be established a month or more before they are mated (Strassmann et al \u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e2004\u003c/span\u003e), presumably due to a scarcity of early males. The unmated replacement queens analysed here are therefore representative of those that would be present on wild nests shortly after queen loss.\u003c/p\u003e \u003cp\u003eTo our knowledge, this study represents the first application of a support vector classification approach to behavior-associated transcriptomic data. Using this approach we identified a large group of genes as differing meaningfully between \u003cem\u003ePolistes\u003c/em\u003e castes\u0026mdash;over 10% of annotated genes, a much larger set than those found using standard analytic approaches in this study and others (e.g. Ferreira et al \u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e2013\u003c/span\u003e; Toth et al \u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e2014\u003c/span\u003e; Patalano et al \u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e2015\u003c/span\u003e; Geffre et al \u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e2017\u003c/span\u003e). Standard approaches using packages such as edgeR (Robinson et al \u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e2010\u003c/span\u003e), DESeq2 (Love et al \u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e2014\u003c/span\u003e), or NOISeq (Tarazona et al \u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e2011\u003c/span\u003e), which typically include information-sharing between genes and relatively strict fold change cutoffs, allow differential expression to be assessed with a high degree of confidence at the level of individual genes. However, because such approaches assess genes individually and do not reduce the dimensionality of the samples\u0026rsquo; transcriptional profiles, they are not well-suited to the tracking of subtle but consistent broad-scale changes. This is especially true for gene expression data that are noisy and heterogenous, or those that are dominated by many genes of small effect (Huang et al \u003cspan citationid=\"CR21\" class=\"CitationRef\"\u003e2018\u003c/span\u003e). The molecular bases of caste in simple social insect societies represent such a class of data.\u003c/p\u003e \u003cp\u003eUsing a machine learning approach, here we have undertaken a detailed analysis of the relationship between gene expression and socially-mediated phenotypic plasticity, revealing broad-scale changes in caste-associated gene expression profiles following a major social disruption. Our results reveal a hitherto unrecognized capacity for large scale disruption to caste-biased gene expression profiles even in the absence of apparent changes in caste phenotype, a disconnect that undermines simplistic models of the relationship between transcriptome and phenotype. Future studies should continue to marry detailed phenotypic and gene expression data in order to assess the prevalence and provenance of such discontinuities.\u003c/p\u003e "},{"header":"Supplementary Methods","content":"\u003ch2\u003eGene expression quantification\u003c/h2\u003e\n\u003cp\u003eBrain tissue was extracted from the heads of individual samples and RNA was extracted using the RNeasy Mini Kit (Qiagen) according to manufacturer\u0026rsquo;s instructions. Library preparation was performed by Novogene Co. followed by sequencing on an Illumina HiSeq 2000 platform with 150-base pair paired-end reads. Reads were filtered with SortMeRNA (Kopylova et al 2012) using default options to remove ribosomal sequences. Trimmomatic (Bolger et al 2014) was then used to perform quality trimming. First, we trimmed adapter sequences and leading and trailing bases with low phred scores (\u0026lt;3). We then used the MAXINFO option with target length 36 and strictness 0.7 to trim low-quality sequences from the remaining reads. Reads were next mapped to 11313 transcripts from the \u003cem\u003eP. dominula\u003c/em\u003e genome annotation 1.0 (Standage et al. 2016) using STAR (Dobin et al 2013) with default options. All 101 samples produced \u0026gt;85% uniquely mapped reads. Reads were assembled into transcripts using StringTie2 (Kovaka et al 2019) before being passed on as raw counts to downstream analyses. Prior to downstream analysis, transcripts were filtered to remove any gene which did not have \u0026gt;20 counts across all samples in at least one of the experimental groups (queens, control workers, and day three and day twelve post-manipulation workers). Following this filtering, 10734/11313 (94.9%) of transcripts remained.\u003c/p\u003e\n\u003ch2\u003eSupport vector machine classification\u003c/h2\u003e\n\u003cp\u003eSupport vector classification was performed in R using the e1071 package (Meyer et al 2017) using the gene expression profiles for all available queens (coded with a value of 1; n=26) and control workers (coded with a value of 0; n=12). Classifiers were assessed via their 3-fold cross-validation error rates; classifiers with lower classification error were considered to be superior. Initially, we tested classifiers using a variety of kernel functions (radial, linear, sigmoid and polynomial) combined with grid searches across a wide range of kernel parameters. A radial kernel with \u0026gamma;=10\u003csup\u003e-6\u003c/sup\u003e and cost parameter C=2\u003csup\u003e5\u003c/sup\u003e was found to produce the lowest error rate of all combinations, so all subsequent classifiers were fit using this kernel and a more focused grid search in a parameter space of\u0026nbsp; 2\u003csup\u003e4\u003c/sup\u003e\u0026lt;C\u0026lt;2\u003csup\u003e6\u003c/sup\u003e and 10\u003csup\u003e-5\u003c/sup\u003e\u0026lt;\u0026gamma;\u0026lt;10\u003csup\u003e-7\u0026nbsp; \u003c/sup\u003ein order to minimise the processing power necessary to perform feature selection. For feature selection, we took the classifier fitted with all genes, and iteratively performed the following process: (1) The 3-fold cross-validation error of the model was calculated twenty times using randomly-assigned bins, and the mean of the resulting errors was recorded as the true validation error of the classifier; (2) The feature weights of all genes were calculated by taking the matrix product of that classifier\u0026rsquo;s coefficients with its support vectors; (3) The gene with the smallest absolute weight in the model was dropped; (4) A new classifier was calculated using the remaining set of genes. This process was repeated until just 100 genes remained, and the optimal support vector classifier was then taken as that for which the cross-validation error reached its minimum.\u003c/p\u003e\n\u003ch2\u003eDifferential expression analysis\u003c/h2\u003e\n\u003cp\u003eDifferential expression analyses were performed in R using the DESeq2 package (Love et al 2014). DESeq2 was run on all groups and contrasts were then calculated for each pair of groups. Differential expression was calculated relative to a baseline fold change of 1.5, i.e. \u003cem\u003ep\u003c/em\u003e-values refer to the probability that absolute change between two groups was greater than 50%. Genes were considered differentially expressed between conditions if \u003cem\u003ep\u003c/em\u003e\u0026lt;0.05 after false discovery rate correction according to the Benjamini-Hochberg procedure.\u003c/p\u003e\n\u003ch2\u003eGene co-expression network analysis\u003c/h2\u003e\n\u003cp\u003eWeighted gene co-expression network analysis was performed in R using the WGCNA package (Langfelder \u0026amp; Horvath 2008). As WGCNA is particularly sensitive to genes with low expression, data were first subjected to a second round of filtering in which genes that had \u0026lt;10 reads in \u0026gt;90% of samples were removed. Consensus gene modules across all samples were then constructed using a soft-threshold power of 9. Initially, 26 gene modules were identified. Modules whose eigengene correlation was \u0026gt;0.75 were subsequently merged, after which 22 consensus modules remained. Finally, the Pearson correlation of each module with each phenotypic trait within each group (queens, control workers and individuals from queen removal nests) was calculated and subjected to Benjamini-Hochberg FDR correction.\u003c/p\u003e\n\u003ch2\u003eGene ontology (GO) enrichment analysis\u003c/h2\u003e\n\u003cp\u003eIn order to perform GO enrichment analysis, we first used OrthoFinder (Emms \u0026amp; Kelly 2019) to identify orthologues for each \u003cem\u003eP. dominula\u003c/em\u003e gene in \u003cem\u003eD. melanogaster\u003c/em\u003e, a model species for which GO annotations are much more complete. GO annotations for each \u003cem\u003eD. melanogaster\u003c/em\u003e gene were acquired from BioMart (Smedley et al 2015) and each \u003cem\u003eP. dominula \u003c/em\u003egene was then assigned GO terms permissively, i.e. a given \u003cem\u003eP. dominula\u003c/em\u003e gene was assigned a GO term if that term appeared as an annotation to any of its orthologues. GO enrichment analysis was then performed in R via the topGO package (Alexa \u0026amp; Rahnenfuhrer 2009) using TopGO\u0026rsquo;s weight01 algorithm and Fisher\u0026rsquo;s exact test to identify GO terms that were significantly overrepresented (\u003cem\u003ep\u003c/em\u003e\u0026lt;0.01) in a focal set of genes against a background consisting of all genes that appeared in the relevant analysis.\u003cbr /\u003e \u003c/p\u003e"},{"header":"Declarations","content":" \u003ch2\u003eAuthor contributions\u003c/h2\u003e \u003cp\u003eBAT, SS and MR conceived of and designed the experiment. BAT and AC carried out the initial experiment. BAT performed behavioural analyses, dissections and RNA extractions, and wrote the initial draft of the manuscript. BAT conducted the analyses with advice from CDRW, MR and SS. BAT and CDRW designed the figures. All authors read and contributed to the manuscript.\u003c/p\u003e \u003ch2\u003eAcknowledgements\u003c/h2\u003e \u003cp\u003eWe would like to thank F. Cappa, I. Pepiciello, L. Dapporto and others at the Dipartimento di Biologia, Universit\u0026agrave; di Firenze for assistance with field work, and E. Favreau and M. Bentley for their assistance with analyses. Special thanks to R. Cervo, without whose generous support this study would not have been possible. This work was supported by the Natural Environment Research Council (grant code NE/L002485/1).\u003c/p\u003e \u003ch2\u003eData availability\u003c/h2\u003e\n\u003cp\u003eSequencing data associated with this paper have been deposited in the NCBI Gene Expression Omnibus (www.ncbi.nlm.nih.gov/geo) under accession number GSE153532 with reviewer access token \u003cem\u003eadctccsahjgjbqt\u003c/em\u003e.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e \u003cspan\u003eAbeel, T., Helleputte, T., Van de Peer, Y., Dupont, P., \u0026amp; Saeys, Y. Robust biomarker identification for cancer diagnosis with ensemble feature selection methods. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cb\u003e26\u003c/b\u003e, 392\u0026ndash;398 (2010).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eAlexa, A., \u0026amp; Rahnenf\u0026uuml;hrer, J. Gene set enrichment analysis with topGO. \u003cem\u003eBioconductor Improv\u003c/em\u003e \u003cb\u003e27\u003c/b\u003e (2009).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eBolger, A. M., Lohse, M., \u0026amp; Usadel, B. Trimmomatic: a flexible trimmer for Illumina sequence data. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cb\u003e30\u003c/b\u003e, 2114\u0026ndash;2120 (2014).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eCasas, L., Saborido-Rey, F., Ryu, T., Michell, C., Ravasi, T., \u0026amp; Irigoien, X. Sex change in clownfish: molecular insights from transcriptome analysis. \u003cem\u003eSci. Rep.\u003c/em\u003e \u003cb\u003e6\u003c/b\u003e, 35461. (2016).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eChakravarty, P., Cozzi, G., Ozgul, A., \u0026amp; Aminian, K. A novel biomechanical approach for animal behaviour recognition using accelerometers. \u003cem\u003eMethods Ecol. Evol.\u003c/em\u003e \u003cb\u003e10\u003c/b\u003e, 802\u0026ndash;814 (2019).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eChandra, V., Fetter-Pruneda, I., Oxley, P. R., Ritger, A. L., McKenzie, S. K., Libbrecht, R., \u0026amp; Kronauer, D. J. Social regulation of insulin signaling and the evolution of eusociality in ants. \u003cem\u003eScience\u003c/em\u003e \u003cb\u003e361\u003c/b\u003e, 398\u0026ndash;402 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eChapman, B. B., Br\u0026ouml;nmark, C., Nilsson, J. \u0026Aring;., \u0026amp; Hansson, L. A. The ecology and evolution of partial migration. \u003cem\u003eOikos\u003c/em\u003e \u003cb\u003e120\u003c/b\u003e, 1764\u0026ndash;1775 (2011).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eCini, A., Meconcelli, S., \u0026amp; Cervo, R. Ovarian indexes as indicators of reproductive investment and egg-laying activity in social insects: a comparison among methods. \u003cem\u003eInsectes Soc.\u003c/em\u003e \u003cb\u003e60\u003c/b\u003e, 393\u0026ndash;402 (2013).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eCullen, D. A., Cease, A. J., Latchininsky, A. V., Ayali, A., Berry, K., Buhl, J., \u0026hellip; Rogers, S. M. From molecules to management: mechanisms and consequences of locust phase polyphenism. \u003cem\u003eAdv. Insect Physiol\u003c/em\u003e \u003cb\u003e53\u003c/b\u003e, 167\u0026ndash;285 (2017).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eDapporto, L., Sledge, F. M., \u0026amp; Turillazzi, S. Dynamics of cuticular chemical profiles of \u003cem\u003ePolistes dominulus\u003c/em\u003e workers in orphaned nests (Hymenoptera, Vespidae). \u003cem\u003eJ. Insect Physiol.\u003c/em\u003e \u003cb\u003e51\u003c/b\u003e, 969\u0026ndash;973 (2005).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eDobin, A., Davis, C. A., Schlesinger, F., Drenkow, J., Zaleski, C., Jha, S., \u0026hellip; Gingeras, T. R. STAR: ultrafast universal RNA-seq aligner. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cb\u003e29\u003c/b\u003e, 15\u0026ndash;21 (2013).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eElo, A. E. (1978). \u003cem\u003eThe rating of chess players, past and present\u003c/em\u003e. New York: Arco Pub.\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eEmms, D. M., \u0026amp; Kelly, S. OrthoFinder: phylogenetic orthology inference for comparative genomics. \u003cem\u003eGenome Biol.\u003c/em\u003e \u003cb\u003e20\u003c/b\u003e, 1\u0026ndash;14 (2019).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eFerreira, P. G., Patalano, S., Chauhan, R., Ffrench-Constant, R., Gabaldon, T., Guigo, R., \u0026amp; Sumner, S. Transcriptome analyses of primitively eusocial wasps reveal novel insights into the evolution of sociality and the origin of alternative phenotypes. \u003cem\u003eGenome Biol.\u003c/em\u003e \u003cb\u003e14\u003c/b\u003e, R20 (2013).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eFowler-Finn, K. D., \u0026amp; Rodr\u0026iacute;guez, R. L. The evolution of experience‐mediated plasticity in mate preferences. \u003cem\u003eJ. Evol. Biol.\u003c/em\u003e \u003cb\u003e25\u003c/b\u003e, 1855\u0026ndash;1863 (2012).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eGeffre, A. C., Liu, R., Manfredini, F., Beani, L., Kathirithamby, J., Grozinger, C. M., \u0026amp; Toth, A. L. Transcriptomics of an extended phenotype: parasite manipulation of wasp social behaviour shifts expression of caste-related genes. \u003cem\u003eProc. R. Soc. B\u003c/em\u003e \u003cb\u003e284\u003c/b\u003e, 20170029 (2017).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eGomulski, L. M., Dimopoulos, G., Xi, Z., Scolari, F., Gabrieli, P., Siciliano, P., \u0026hellip; Gasperi, G. Transcriptome profiling of sexual maturation and mating in the Mediterranean fruit fly, \u003cem\u003eCeratitis capitata\u003c/em\u003e. \u003cem\u003ePLoS One\u003c/em\u003e \u003cb\u003e7\u003c/b\u003e, (2012).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eHadfield, J. D., Nutall, A., Osorio, D., \u0026amp; Owens, I. P. F. Testing the phenotypic gambit: phenotypic, genetic and environmental correlations of colour. \u003cem\u003eJ. Evol. Biol.\u003c/em\u003e \u003cb\u003e20\u003c/b\u003e, 549\u0026ndash;557 (2007).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eHeyes, C. Blackboxing: social learning strategies and cultural evolution. \u003cem\u003ePhilos. Trans. R. Soc. B\u003c/em\u003e \u003cb\u003e371\u003c/b\u003e, 20150369 (2016).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eHoriguchi, R., Nozu, R., Hirai, T., Kobayashi, Y., Nagahama, Y., \u0026amp; Nakamura, M. Characterization of gonadal soma-derived factor expression during sex change in the protogynous wrasse, \u003cem\u003eHalichoeres trimaculatus\u003c/em\u003e. \u003cem\u003eDev. Dynam.\u003c/em\u003e \u003cb\u003e242\u003c/b\u003e, 388\u0026ndash;399 (2013).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eHuang, S., Cai, N., Pacheco, P. P., Narrandes, S., Wang, Y., \u0026amp; Xu, W. Applications of support vector machine (SVM) learning in cancer genomics. \u003cem\u003eCancer Genom. Proteom.\u003c/em\u003e \u003cb\u003e15\u003c/b\u003e, 41\u0026ndash;51 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eJandt, J. M., Tibbetts, E. A., \u0026amp; Toth, A. L. \u003cem\u003ePolistes\u003c/em\u003e paper wasps: a model genus for the study of social dominance hierarchies. \u003cem\u003eInsectes Soc. 61\u003c/em\u003e, 11\u0026ndash;27 (2014).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eKopylova, E., No\u0026eacute;, L., \u0026amp; Touzet, H. SortMeRNA: fast and accurate filtering of ribosomal RNAs in metatranscriptomic data. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cb\u003e28\u003c/b\u003e, 3211\u0026ndash;3217 (2012).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eKovaka, S., Zimin, A. V., Pertea, G. M., Razaghi, R., Salzberg, S. L., \u0026amp; Pertea, M. Transcriptome assembly from long-read RNA-seq alignments with StringTie2. \u003cem\u003eGenome Biol.\u003c/em\u003e \u003cb\u003e20\u003c/b\u003e, 1\u0026ndash;13 (2019).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eLangfelder, P., \u0026amp; Horvath, S. WGCNA: an R package for weighted correlation network analysis. \u003cem\u003eBMC Bioinformatics\u003c/em\u003e \u003cb\u003e9\u003c/b\u003e, 559 (2008).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eLibbrecht, R., Oxley, P. R., \u0026amp; Kronauer, D. J. Clonal raider ant brain transcriptomics identifies candidate molecular mechanisms for reproductive division of labor. \u003cem\u003eBMC Biol.\u003c/em\u003e \u003cb\u003e16\u003c/b\u003e, 89 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eLo, N., Simpson, S. J., \u0026amp; Sword, G. A. Epigenetics and developmental plasticity in orthopteroid insects. \u003cem\u003eCurr. Opin. Insect. Sci.\u003c/em\u003e \u003cb\u003e25\u003c/b\u003e, 25\u0026ndash;34 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eLove, M., Anders, S., \u0026amp; Huber, W. Differential analysis of count data\u0026ndash;the DESeq2 package. \u003cem\u003eGenome Biol.\u003c/em\u003e \u003cb\u003e15\u003c/b\u003e, 10\u0026ndash;1186 (2014).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eManfredini, F., Brown, M. J., \u0026amp; Toth, A. L. Candidate genes for cooperation and aggression in the social wasp \u003cem\u003ePolistes dominula\u003c/em\u003e. \u003cem\u003eJ. Comp. Physiol. A\u003c/em\u003e \u003cb\u003e204\u003c/b\u003e, 449\u0026ndash;463 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eMatson, C. K., Murphy, M. W., Sarver, A. L., Griswold, M. D., Bardwell, V. J., \u0026amp; Zarkower, D. DMRT1 prevents female reprogramming in the postnatal mammalian testis. \u003cem\u003eNature\u003c/em\u003e \u003cb\u003e476\u003c/b\u003e, 101\u0026ndash;104 (2011).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eMeyer, D., Dimitriadou, E., Hornik, K., Weingessel, A. \u0026amp; Leisch, F. e1071: Misc Functions of the Department of Statistics, Probability Theory Group (Formerly: E1071), TU Wien. R package version 1.6-8. \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://CRAN.R-project.org/package=e1071\u003c/span\u003e\u003c/span\u003e. (2017).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eMiyano, S. Colony development, worker behavior and male production in orphan colonies of a Japanese paper wasp, \u003cem\u003ePolistes chinensis antennalis\u003c/em\u003e P\u0026eacute;rez. \u003cem\u003ePopul. Ecol.\u003c/em\u003e \u003cb\u003e68\u003c/b\u003e, 347\u0026ndash;361 (1986).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eMiyano, S. Worker reproduction and related behavior in orphan colonies of a Japanese paper wasp, \u003cem\u003ePolistes jadwigae\u003c/em\u003e (Hymenoptera, Vespidae). \u003cem\u003eJ. Ethol.\u003c/em\u003e \u003cb\u003e9\u003c/b\u003e, 135\u0026ndash;146 (1991).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eMonnin, T., Cini, A., Lecat, V., F\u0026eacute;d\u0026eacute;rici, P., \u0026amp; Doums, C. No actual conflict over colony inheritance despite high potential conflict in the social wasp \u003cem\u003ePolistes dominulus\u003c/em\u003e. \u003cem\u003eProc. R. Soc. B\u003c/em\u003e \u003cb\u003e276\u003c/b\u003e, 1593\u0026ndash;1601 (2009).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eNeumann, C., Duboscq, J., Dubuc, C., Ginting, A., Irwan, A. M., Agil, M., \u0026hellip; Engelhardt, A. Assessing dominance hierarchies: validation and advantages of progressive evaluation with Elo-rating. \u003cem\u003eAnim. Behav.\u003c/em\u003e \u003cb\u003e82\u003c/b\u003e, 911\u0026ndash;921 (2011).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003ePardi, L. Dominance order in \u003cem\u003ePolistes\u003c/em\u003e wasps. \u003cem\u003ePhysiol. Zool.\u003c/em\u003e \u003cb\u003e21\u003c/b\u003e, 1\u0026ndash;13 (1948).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003ePatalano, S., Vlasova, A., Wyatt, C., Ewels, P., Camara, F., Ferreira, P. G., \u0026hellip; Sumner, S. Molecular signatures of plastic phenotypes in two eusocial insect species with simple societies. \u003cem\u003eProc. Natl. Acad. Sci. U.S.A.\u003c/em\u003e \u003cb\u003e112\u003c/b\u003e, 13970\u0026ndash;13975 (2015).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eRast, W., Kimmig, S. E., Giese, L., \u0026amp; Berger, A. Machine learning goes wild: Using data from captive individuals to infer wildlife behaviours. \u003cem\u003ePLoS One\u003c/em\u003e \u003cb\u003e15\u003c/b\u003e, e0227317 (2020).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eR\u0026eacute;ale, D., Dingemanse, N. J., Kazem, A. J., \u0026amp; Wright, J. Evolutionary and ecological approaches to the study of personality. \u003cem\u003ePhilos. Trans. R. Soc. B.\u003c/em\u003e \u003cb\u003e365\u003c/b\u003e, 3947\u0026ndash;3946 (2010).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eRehan, S. M., Glastad, K. M., Steffen, M. A., Fay, C. R., Hunt, B. G., \u0026amp; Toth, A. L. Conserved genes underlie phenotypic plasticity in an incipiently social bee. \u003cem\u003eGenome Biol. Evol.\u003c/em\u003e \u003cb\u003e10\u003c/b\u003e, 2749\u0026ndash;2758 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eRittschof, C. C., \u0026amp; Robinson, G. E. Genomics: moving behavioural ecology beyond the phenotypic gambit. \u003cem\u003eAnim. Behav.\u003c/em\u003e \u003cb\u003e92\u003c/b\u003e, 263\u0026ndash;270 (2014).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eRobinson, M. D., McCarthy, D. J., \u0026amp; Smyth, G. K. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. \u003cem\u003eBioinformatics\u003c/em\u003e \u003cb\u003e26\u003c/b\u003e, 139\u0026ndash;140 (2010).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eRubin, H. The phenotypic gambit: selective pressures and ESS methodology in evolutionary game theory. \u003cem\u003eBiol. Philos.\u003c/em\u003e \u003cb\u003e31\u003c/b\u003e, 551\u0026ndash;569 (2016).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eSchwander, T., Lo, N., Beekman, M., Oldroyd, B. P., \u0026amp; Keller, L. Nature versus nurture in social insect caste differentiation. \u003cem\u003eTrends Ecol. Evol.\u003c/em\u003e \u003cb\u003e25\u003c/b\u003e, 275\u0026ndash;282 (2010).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eSegal, N. H., Pavlidis, P., Antonescu, C. R., Maki, R. G., Noble, W. S., DeSantis, D., \u0026hellip; Cordon-Cardo, C. Classification and subtype prediction of adult soft tissue sarcoma by functional genomics. \u003cem\u003eAm. J. Pathol.\u003c/em\u003e \u003cb\u003e163\u003c/b\u003e, 691\u0026ndash;700 (2003).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eShell, W. A., \u0026amp; Rehan, S. M. Behavioral and genetic mechanisms of social evolution: insights from incipiently and facultatively social bees. \u003cem\u003eApidologie\u003c/em\u003e \u003cb\u003e49\u003c/b\u003e, 13\u0026ndash;30 (2018).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eSimola, D. F., Graham, R. J., Brady, C. M., Enzmann, B. L., Desplan, C., Ray, A., \u0026hellip; Berger, S. L. Epigenetic (re) programming of caste-specific behavior in the ant \u003cem\u003eCamponotus floridanus\u003c/em\u003e. \u003cem\u003eScience\u003c/em\u003e \u003cb\u003e351\u003c/b\u003e, aac6633 (2016).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eSmedley, D., Haider, S., Durinck, S., Pandini, L., Provero, P., Allen, J., \u0026hellip; Bardou, P. The BioMart community portal: an innovative alternative to large, centralized data repositories. \u003cem\u003eNucleic Acids Res.\u003c/em\u003e \u003cb\u003e43\u003c/b\u003e, W589-W598 (2015).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eStandage, D. S., Berens, A. J., Glastad, K. M., Severin, A. J., Brendel, V. P., \u0026amp; Toth, A. L. Genome, transcriptome and methylome sequencing of a primitively eusocial wasp reveal a greatly reduced DNA methylation system in a social insect. \u003cem\u003eMol. Ecol.\u003c/em\u003e \u003cb\u003e25\u003c/b\u003e, 1769\u0026ndash;1784 (2016).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eStarks, P. T., Turillazzi, S., \u0026amp; West-Eberhard, M. J. \u003cem\u003ePolistes\u003c/em\u003e paper wasps: emergence of a model genus. In \u003cem\u003eAnn. Zool. Fenn.\u003c/em\u003e \u003cb\u003e43\u003c/b\u003e, 385\u0026ndash;386 (2006).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eStrassmann, J. E., Fortunato, A., Cervo, R., Turillazzi, S., Damon, J. M., \u0026amp; Queller, D. C. The cost of queen loss in the social wasp \u003cem\u003ePolistes dominulus\u003c/em\u003e (Hymenoptera: Vespidae). \u003cem\u003eJ. Kansas Entomol. Soc.\u003c/em\u003e \u003cb\u003e77\u003c/b\u003e, 343\u0026ndash;355 (2004).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eTarazona, S., Garc\u0026iacute;a, F., Ferrer, A., Dopazo, J., \u0026amp; Conesa, A. NOIseq: a RNA-seq differential expression method robust for sequencing depth biases. \u003cem\u003eEMBnet.journal\u003c/em\u003e \u003cb\u003e17\u003c/b\u003e, 18\u0026ndash;19. (2011).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eTaylor, B. A., Cini, A., Cervo, R., Reuter, M., \u0026amp; Sumner, S. Queen succession conflict in the paper wasp \u003cem\u003ePolistes dominula\u003c/em\u003e is mitigated by age-based convention. \u003cem\u003eBehav. Ecol.\u003c/em\u003e, araa045 (2020).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eTibbetts, E. A., \u0026amp; Huang, Z. Y. The challenge hypothesis in an insect: juvenile hormone increases during reproductive conflict following queen loss in \u003cem\u003ePolistes\u003c/em\u003e wasps. \u003cem\u003eAm. Nat.\u003c/em\u003e \u003cb\u003e176\u003c/b\u003e, 123\u0026ndash;130 (2010).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eTodd, E. V., Ortega-Recalde, O., Liu, H., Lamm, M. S., Rutherford, K. M., Cross, H., \u0026hellip; Godwin, J. R. Stress, novel sex genes, and epigenetic reprogramming orchestrate socially controlled sex change. \u003cem\u003eSci. Adv.\u003c/em\u003e \u003cb\u003e5\u003c/b\u003e, eaaw7006 (2019).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eToth, A. L., Tooker, J. F., Radhakrishnan, S., Minard, R., Henshaw, M. T., \u0026amp; Grozinger, C. M. Shared genes related to aggression, rather than chemical communication, are associated with reproductive dominance in paper wasps (\u003cem\u003ePolistes metricus\u003c/em\u003e). \u003cem\u003eBMC Genom.\u003c/em\u003e \u003cb\u003e15\u003c/b\u003e, 75 (2014).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eTsuji, K., \u0026amp; Tsuji, N. Why is dominance hierarchy age-related in social insects? The relative longevity hypothesis. \u003cem\u003eBehav. Ecol. Sociobiol.\u003c/em\u003e \u003cb\u003e58\u003c/b\u003e, 517\u0026ndash;526 (2005).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eWest-Eberhard, M. J. The Social Biology of Polistine Wasps. \u003cem\u003eMisc. Pub. Mus. Zool. Univ. Mich.\u003c/em\u003e \u003cb\u003e140\u003c/b\u003e, 1\u0026ndash;101 (1969).\u003c/span\u003e \u003c/li\u003e \u003cli\u003e \u003cspan\u003eZhou, S., Mackay, T. F., \u0026amp; Anholt, R. R. Transcriptional and epigenetic responses to mating and aging in \u003cem\u003eDrosophila melanogaster\u003c/em\u003e. \u003cem\u003eBMC Genom.\u003c/em\u003e \u003cb\u003e15\u003c/b\u003e, 927 (2014).\u003c/span\u003e \u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":false,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"nature-portfolio","isNatureJournal":true,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"","title":"Nature Portfolio","twitterHandle":"","acdcEnabled":false,"dfaEnabled":false,"editorialSystem":"ejp","reportingPortfolio":"","inReviewEnabled":true,"inReviewRevisionsEnabled":false},"keywords":"Phenotypic plasticity, machine learning, Polistes dominula","lastPublishedDoi":"10.21203/rs.3.rs-43815/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-43815/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003ePhenotypic plasticity, the ability to produce multiple phenotypes from a single genotype, represents an excellent model with which to examine the relationship between gene expression and phenotypes. Despite this, analyses of the molecular bases of plasticity have been limited by the challenges of linking individual phenotypes with individual-level gene expression profiles, especially in the case of complex social phenotypes. Here, we tackle this challenge by analysing the individual-level gene expression profiles of \u003cem\u003ePolistes dominula\u003c/em\u003e paper wasps following the loss of a queen, a perturbation that induces some individuals to undergo a significant phenotypic shift and become replacement reproductives. Using a machine learning approach, we find a strong response of caste-associated gene expression to queen loss, wherein individuals\u0026rsquo; expression profiles become intermediate between queen and worker states. Importantly, this change occurs even in individuals that appear phenotypically unaffected. Part of this response is explained by individual attributes, most prominently age. These results demonstrate that large changes in gene expression may occur in the absence of detectable phenotypic changes, resulting here in a socially mediated de-differentiation of individuals at the transcriptomic but not the phenotypic level. Our findings also highlight the complexity of the relationship between gene expression and phenotype, where transcriptomes are neither a direct reflection of the genotype nor a proxy for the molecular underpinnings of the external phenotype.\u003c/p\u003e","manuscriptTitle":"The molecular basis of socially-mediated phenotypic plasticity in a eusocial paper wasp","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2020-07-27 18:43:59","doi":"10.21203/rs.3.rs-43815/v1","editorialEvents":[],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"nature-communications","isNatureJournal":true,"hasQc":false,"allowDirectSubmit":false,"externalIdentity":"NCOMMS","sideBox":"Learn more about [Nature Communications](http://www.nature.com/ncomms/)","snPcode":"","submissionUrl":"https://mts-ncomms.nature.com/","title":"Nature Communications","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"ejp","reportingPortfolio":"Nature Communications","inReviewEnabled":true,"inReviewRevisionsEnabled":false}}],"origin":"","ownerIdentity":"bad3bd00-7ba8-4114-991a-482d20072901","owner":[],"postedDate":"July 27th, 2020","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"published-in-journal","subjectAreas":[{"id":207557,"name":"Entomology"},{"id":207558,"name":"Animal Behavior"},{"id":207559,"name":"Bioinformatics"}],"tags":[],"updatedAt":"2021-08-20T02:42:17+00:00","versionOfRecord":{"articleIdentity":"rs-43815","link":"https://doi.org/10.1038/s41467-021-21095-6","journal":{"identity":"nature-communications","isVorOnly":false,"title":"Nature Communications"},"publishedOn":"2021-02-03 02:42:17","publishedOnDateReadable":"February 3rd, 2021"},"versionCreatedAt":"2020-07-27 18:43:59","video":"","vorDoi":"10.1038/s41467-021-21095-6","vorDoiUrl":"https://doi.org/10.1038/s41467-021-21095-6","workflowStages":[]},"version":"v1","identity":"rs-43815","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-43815","identity":"rs-43815","version":["v1"]},"buildId":"7rjqhiLT3MXkJMwkYKINL","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}