Native and introduced Argentine ant populations are characterised by distinct transcriptomic signatures associated with behaviour and immunity

populations Abstract Biological invasions can be influenced by trait variation in the invader, such as behavioural traits and ecological factors, such as variation in pathogen pressure. High-throughput nucleotide sequencing has increased our capacity to investigate the genomic basis of the functional changes associated with biological invasions. Here, we used RNA-sequencing in Argentina and California, Australia and New Zealand to investigate if native and introduced Argentine ant populations were characterised by distinct transcriptomic signatures. We focused our analysis on viral pressure and immunity, as well as genes associated with biogenic amines known to modulate key behaviour in social insects. Using a combination of differential expression analysis, gene co-expression network analysis and candidate gene approach, we show that native and introduced populations have distinct transcriptomic signatures. Genes associated with biogenic amines were overall up-regulated in the native range compared to introduced populations. Although we found no significant variation in overall viral loads amongst regions for viruses known to infect Argentine ants, viral diversity was lower in most of the introduced range which was interestingly associated with down-regulation of the RNAi immune pathway, primarily directed against viruses. Altogether, our data show that Argentine ant populations exhibit range-specific transcriptomic signatures, perhaps reflecting regional adaptations that may contribute to the ecological success of introduced populations.


Introduction
Exotic species are commonly transported around the world as inadvertent stowaways in cargo and can sometime become invasive and pose great threats to biodiversity, agriculture and other human activities (Hulme 2009).Through population bottlenecks or novel evolutionary forces, both during the introduction process and in the introduced range, introduced species can experience rapid changes and sometimes enhanced ecological success (Sakai et al. 2001;Suarez and Tsutsui 2008).Some traits and ecological factors, such as behavioural variation and pathogen pressure, are key determinants of invasion success and, therefore, quantifying variation in these may be particularly relevant for understanding biological invasions.High-throughput sequencing has increased our capacity to investigate the genomic basis of the functional changes associated with biological invasions (Rius et al. 2015).Here, we used RNA-sequencing data to investigate possible functional variations between native and introduced populations of the Argentine ant (Linepithema humile), a globally successful pest.We focused our analysis on genes associated with biogenic amines, which modulate key behaviour in social insects (Kamhi et al. 2017), as well as immunity and associated virus diversity.
Some behavioural traits are regarded as major drivers of biological invasions and have been specifically suggested in the context of ant invasions (Holway and Suarez 1999;Phillips and Suarez 2012;Silverman and Buczkowski 2016).For instance, low intraspecific aggression is thought to be a key driver of the ecological success of the Argentine ant in its introduced range (Holway et al. 1998).Other behavioural traits, such as increased foraging activity and high interspecific aggression, are also likely to influence the Argentine ant's invasiveness (Silverman and Buczkowski 2016).Neural pathways associated with hallmark behaviour of social insects, such as foraging and aggression, have been shown to be modulated by compounds, including the biogenic amines octopamine (OA), serotonin (5-HT) and dopamine (DA) (Kamhi et al. 2017).Interestingly, such behavioural traits are often correlated and can be considered as components of behavioural syndromes (Jandt et al. 2014).Biogenic amines can function as neurotransmitters (neuron-to-neuron communication), but can also function as neuromodulators and neurohormones, in which case they target larger regions of the brain or the body, respectively, and underlie specific behaviour (Hoyle 1985;Libersat and Pflueger 2004).In the context of biological invasions, variation in biogenic amine pathways could be associated with the regulation of behavioural syndromes that contribute to the success of introduced species.In social insects, expression of genes associated with biogenic amines, has been suggested to be associated with variation in social behaviour, including foraging and aggression (Liang et al. 2012(Liang et al. , 2014;;Kamhi et al. 2017;Friedman et al. 2018).Some of these genes have been shown to exhibit genetic variation associated with behavioural syndromes (e.g. in birds: Fidler et al. 2007;Korsten et al. 2010) and can exhibit range-specific polymorphism (Mueller et al. 2013;Holtmann et al. 2016).Behavioural differences between the native and introduced Argentine ant populations might be restricted to specific behaviour in a specific context (Felden et al. 2018).To gain further insight into variation in the molecular basis of the Argentine ant's behaviour across its native and introduced ranges, we measured expression of genes related to the main biogenic amine pathways.We hypothesised that the molecular basis of behaviour, central to the fitness of ant societies, may exhibit rangespecific expression profiles.
Pathogens often play an important role in biological invasions (Tompkins et al. 2011;Dunn et al. 2012;White and Perkins 2012).Social insects possess several characteristics that likely make them vulnerable to pathogens, including high densities of individuals within nests and social groups comprised of related individuals (Schmid-Hempel 1998).Introduced species, including invasive ants, can exhibit boom-and-bust population dynamics and pathogens may play a part in such dynamics (Simberloff and Gibbons 2004;Lester and Gruber 2016).The Argentine ant is known to harbour a number of microbes and invasion pathways have been shown to be associated with changes in microbial communities of invasive ants, including both losses and gains of pathogens and endosymbionts (Yang et al. 2010;Sébastien et al. 2015;Gruber et al. 2017;Lester et al. 2017;Viljakainen et al. 2018).Either loss or acquisition of pathogens along an introduction pathway could directly affect an invader's fitness and influence population dynamics and invasion success.Furthermore, under the Evolution of Increased Competitive Ability hypotheses (EICA), diverting resources from immunity or allocating resources to less costly immune pathways in response to relaxed pathogen pressure is thought to promote invasion success (Keane and Crawley 2002;Torchin et al. 2003;Lee and Klasing 2004).Homology-based analysis has revealed many immune genes that appear conserved across model organisms and Argentine ants (Smith et al. 2011).We measured variation in the expression of genes associated with the JaK/STAT, RNAi, Toll and Imd immune signalling pathways across the native and introduced ranges.To clarify the landscape of viral pressure in the native and invaded ranges, we also measured viral diversity and overall viral loads in our samples.We hypothesised that, if Argentine ant invasion were facilitated by variation in viral pressure, introduced populations may exhibit differential regulation of immune genes, thereby promoting invasiveness both directly and through beneficial resource re-allocation.
We examined possible functional adaptations underlying the success of a globally invasive pest, using RNA-sequencing across the native and introduced ranges of the Argentine ant.Are introduced populations characterised by differences in gene expression that could underpin behavioural variation?Are introduced populations characterised by a release from pathogens from the native range and does that translate into lower immune gene expression?To answer these questions, we investigated 1) variation in expression of genes associated with biogenic amines and 2) immunity, as well as 3) viral diversity in the native range (Argentina) and the introduced range in California, Europe, Australia and New Zealand.

Sampling and RNA library preparation
We used worker ants collected in Argentina, California, Australia and New Zealand from colonies maintained in standardised conditions for 20 days prior to sampling as described in Felden et al. (2018).Ants were collected at four different sites in each region, except in Europe that only included two distinct sites (Suppl.material 1: Table S1).Briefly, colonies were maintained in experimental colonies comprising 1,200 workers and four queens, fed daily with one mealworm (Tenebrio molitor, ≈ 150 mg) and 1 ml of 20% sucrose solution (or 1% sucrose for Low sugar colonies, i.e. samples AR OTA-2 , CA RCH-2 , AU STA-1 , NZ HAS-3 and NZ HAS-4 ).For the purpose of virus presence and load analysis only (i.e.not gene expression), we also included samples from colonies that were fed with one mealworm and octopamine in 1% sucrose solution (i.e.samples AR OTA-3 , CA RCH-3 , AU STA-2 , AU STA-3 ; see Felden et al. 2018).Foraging workers were sampled and stored in ice-cold RNALater (Invitrogen, USA).Tubes of ants in RNALater were kept at 6 °C for 24 h, and at -20 °C in the country of origin until shipped to New Zealand where they were stored at -80 °C for up to six months until RNA was extracted.
Ant heads with antennae were separated from bodies in RNALater under a stereomicroscope and total RNA was extracted from 20 pooled heads and antennae of workers from the same colony.Samples were then briefly washed with ice-cold phosphate-buffered saline (PBS) to remove RNALater that can affect extraction quality.RNA was extracted using an in-column Trizol-based purification kit, using the manufacturer's recommended methods (Direct-Zol Microprep, Zymo Research, USA).RNA integrity was confirmed and quantified with an RNA 6000 Nano chip on the Agilent 2100 Bioanalyzer (Agilent Technologies Co. Ltd., Diegem, Belgium), according to the manufacturer's instructions.Extracted RNA was stored in RNAStable (Biomatrica Inc., San Diego, USA) and sent to BGI (Shenzen, China) for Illumina Hi-Seq sequencing.Overall, 30 head/antenna libraries were sequenced in the native and introduced ranges, including six replicates from Argentina (four sites), seven from California (four sites), three in Europe (two sites), seven from Australia (four sites) and seven from New Zealand (three sites) (see details in Suppl.material 1: Table S1).Four libraries that were experimentally treated with octopamine (OA) for another study were discarded from the gene expression analysis and only used in the virus analysis, as experimental treatment likely affected gene expression.Samples were sequenced as 150 base paired-end barcoded mRNA TruSeq libraries, aiming to generate 4 GB of data per sample.Pre-processing at BGI included the removal of reads with more than 10% values missing, reads with more than 10% of bases with quality scores Q < 20 and removal of adapters.

Read alignment and Argentine ant transcript quantification
Computationally demanding analyses were performed on Victoria University of Wellington Science Faculty's High Performance Computing Facility.Clean paired-end reads were aligned to the Argentine ant reference genome (Smith et al. 2011) with HISAT 2.1.0using default parameters (Kim et al. 2015) to produce sample-specific BAM files.Overall, the average read alignment rate to the Argentine ant genome was 89.66%.StringTie 1.3.4(Pertea et al. 2015) was used to generate GTF files from the BAM files, using the -e argument to restrict the assembly to transcripts matching the reference (Lhum_UMD_V04, accession number GCF_000217595.1).We generated a raw transcript counts matrix at the gene level using StringTie's authors' prepDE.pyscript.

Differential gene expression analysis
To detect the most differentially expressed genes between the native and the introduced ranges, we followed part of the guidelines outlined in Law et al. (2016).Briefly, we imported the raw counts matrix into R as a DGEList object using the edgeR 3.22.3package (Robinson et al. 2010), together with one of the GTF files created by StringTie as a source of information for each gene.To filter the raw counts matrix for low expressed transcripts, we first computed counts per million (CPM), discarded transcripts expressed in less than three samples at more than 1 CPM and then applied a TMM normalisation on the filtered raw count data.We then used the limma 3.36.5 (Ritchie et al. 2015) in R to perform a series of four pairwise differential gene expression analyses between Argentina and each region of the introduced range included in the dataset (i.e.California, Europe, Australia, New Zealand).We transformed TMM-normalised raw counts filtered for low expressed transcripts into log2-counts per million (log2CPM) with the limma/voom function and fitted a linear model for each gene on this post-filtering TMMnormalised log2CPM matrix.We selected the most differentially expressed genes within each pairwise comparison using a FDR cut-off of 0.05 and an absolute fold change > 1.1.To identify genes consistently differentially expressed between the native and the introduced range, we selected the overlap of the differentially expressed gene in all four pairwise comparisons, which yielded a set of 130 genes.In order to gain a broader functional assessment of the lists of differentially expressed genes between the native and introduced ranges, we undertook gene ontology (GO) enrichment analysis.First, we used BLASTx to determine the closest matches of all Argentine ant reference genes (assembly Lhum_ UMD_V04, accession number GCF_000217595.1) against the Swissprot database (The Uniprot Consortium 2017; downloaded 28/02/2019).Then, we imported the BLASTx output as a XML file into BLAST2GO and proceeded to the mapping and annotation step in order to run Gene Set Enrichment Analyses (GSEA) restricted to Biological Functions.The list of differentially expressed genes was ranked based on fold-change.

Weighted Gene Co-expression Network Analysis (WGCNA)
WGCNA (in WGCNA 1.64-1, Langfelder and Horvath 2008) allows the detection of modules of co-expressed genes that can be correlated to factors such as pheno-typic traits and is an alternative approach to detect genes of interest to that of a classical differential gene expression analysis as described above.To examine transcriptome-wide expression patterns associated with range (native or introduced), we analysed the post-filtering TMM-normalised log2CPM matrix produced prior to the differential gene expression analysis using WGCNA.We followed WGCNA authors' guidelines to further detect low expressed transcripts and outlying samples.We then used the scale-free topology criterion to select a soft thresholding power to build the gene co-expression network (Langfelder and Horvath 2008; Suppl.material 1: Figure S2).We used range as a numerical variable (0 for native populations and 1 for introduced populations) to investigate the correlation between the eigengene of each module and range.We investigated GO enrichment of the module that was the most correlated with range as outlined above in the DGE analysis section, except that genes were ranked based on module membership.Module membership indicates how the focal gene expression is correlated with the module eigengene expression profile (Langfelder and Horvath 2008).Highly connected genes that exhibit high module membership are likely to have particular biological relevance.

Candidate gene approach
To determine differences in gene expression associated with biogenic amines and immunity, we compiled a list of genes of interest, based on existing annotations of the Argentine ant genome (Lhum_UMD_V04 from Smith et al. 2011, listed in Suppl.material 1: Table S3, Table S4, respectively).First, we computed log-centred FPKM counts for genes of interest using the rpkm function in edgeR (Robinson et al. 2010), using the post-filtering TMM-normalised counts generated from the DGEList object described in the differential gene expression analysis section (i.e.prior to log2-CPM transformation).For genes belonging to immune signalling pathways, we computed the average expression per region (i.e.Argentina or one of the introduced regions) per gene and tested for differences in overall immune pathway expression.For immune pathways that contained 10 annotated genes or more, we analysed overall immune pathway gene expression, using data for all genes associated with each pathway of interest using linear models, with region as an independent variable.We tested for post-hoc pairwise differences amongst groups using the multcomp/glht R function (Hothorn et al. 2016).For immune pathways that comprised less than 10 annotated genes (i.e.JaK/STAT and JNK), we used a Kruskal-Wallis test, followed by post-hoc multiple comparisons if significant (Giraudoux et al. 2018).The genes of interest, associated with biogenic amines, were in limited number (i.e. 3 to 4 genes per biogenic amine pathway), hence we restricted the statistical analysis at the gene level and compared variation associated with range by pooling the data from the introduced regions.

Virome analysis
Reads that did not align to the Argentine ant genome were assembled de novo using Trinity 2.3.2 (Haas et al. 2013) using default parameters.We quantified transcript expression within Trinity 2.3.2 using eXpress (Roberts and Pachter 2012), yielding a TMM-normalised TPM matrix at the gene level.We then aligned the assembled transcripts to various reference databases to screen for and quantify viral transcripts.We used BLAST 2.2.3 (Camacho et al. 2009) to run BLASTx searches of the NCBI viral protein database (downloaded on 19/11/2018 from ftp://ftp.ncbi.nlm.nih.gov/refseq/release/viral/), as well as of a more recently published Argentine ant novel viruses database (Viljakainen et al. 2018), using a e-value cut-off of 10 -5 .We discarded transcripts that were less than 100 bases long and queries that returned < 90% identity with the target sequence.From the filtered BLAST outputs, we selected a single best hit per assembled transcript based on the highest bitscore.To compute viral loads, we first summed TPMs for all genes belonging to each identified virus (Suppl.material 1: File S1).Then, we used the Argentine ant library size, computed in the differential expression analysis section, to normalise viral loads to host tissue amongst samples.
Virus presence was verified via RT-PCR.Briefly, RNA was extracted from 10 whole ants per site (i.e.four extractions per region) following a similar Direct-Zol protocol as described above.Ants from Europe were sampled at only two different locations; hence these samples included two replicated extractions per site.A total of 250 ng of RNA from each extraction was pooled with respect to region so that each regional sample contained 1 µg of RNA.We prepared cDNA libraries using qScript cDNA SuperMix (Quanta Biosciences), using the manufacturer's instructions.Target-specific PCR conditions are given in Suppl.material 1: Table S6.

Overall transcriptomic profile variation and most differentially expressed genes between the native and introduced ranges
Overall gene expression profiles clustered with region and range (Figure 1).We found highly variable numbers of differentially expressed genes between Argentina and introduced regions (Suppl.material 1: Figure S1, Table S2).For further analysis of rangespecific transcriptomic signatures, we considered the 130 genes that were consistently differentially expressed between Argentina and all of the introduced populations (genes coloured in red in Suppl.material 1: Figure S1, File S1).
Genes that were consistently differentially expressed between the native and introduced ranges tended to be up-regulated in the introduced range (118/130 genes; Figure 2, Suppl.material 1: File S1).Gene set enrichment analysis indicated that this set of differentially expressed genes was only significantly enriched for  S1).al 1: Figure S3).One of these modules was significantly (p < 0.001) and strongly correlated (r = 0.96) with range (Suppl.material 1: Figure S4).A total of 103 genes within the module were also present in the set of differentially expressed genes between the native and introduced ranges (Suppl.material 1: File S1).Gene set enrichment analysis of the module, associated with range, showed significant enrichment for redundant biological process terms related to immune signalling pathways (e.g.positive regulation of response to biotic stimulus, immune response-activating signal transduction, regulation of immune effector process with up to 36 genes per enriched term, FDR < 0.25, Suppl.material 1: File S3) and phototransduction (e.g.sensory perception of light stimulus, visual perception with up to 41 genes per enriched term, FDR < 0.25, Suppl.material 1: File S3).Thirteen genes overlapped between our list of candidate immune genes (Suppl.material 1: Table S4) and the coexpressed genes in the module, including five Toll-like transcripts (LOC105668121, LOC105668599, LOC105674013, LOC105674712 and LOC105678784), four beta-1,3-glucan-binding protein-like transcripts (LOC105673881, LOC105674418, LOC105674426 and LOC105674427), a PGRP-LC transcript (LOC105678063), defensin-2 (LOC105975717), serine protease nudel-like (LOC105670067) and a Toll/ Imd-associated uncharacterised protein transcript (LOC105672003) (The Uniprot Consortium 2017).Five of these transcripts (i.e.Toll-like protein LOC105668121, beta-1,3-glucan-binding protein-like transcripts LOC105674418, LOC105674426, LOC10567442, and uncharacterised LOC105672003) were amongst the 20% of most connected transcripts within the module (Suppl.material 1: Figure S5).

Candidate gene approach: genes associated with biogenic amines
Amongst the 14 genes selected for their association with biogenic amine pathways in Argentine ants, 10 serotoninergic, dopaminergic, octopaminergic and tyraminergic receptors exhibited expression levels that were significantly higher in the native range compared to the introduced range (Figure 3).A transcript coding for tyramine beta-hydroxylase (TyrBH) was also significantly up-regulated in Argentina and only the octopaminergic receptor OctR1 and two serotoninergic receptors (5-HTR-2A and 5-HTR-2C) expression levels were not significantly different between the two ranges.

Candidate gene approach: immune pathways
An immune pathway-specific analysis of gene expression revealed that genes associated with the RNAi pathway were consistently down-regulated in the introduced range compared to Argentina (p < 0.001 in all introduced regions; Figure 4a, Suppl.material 1: Table S5).The JaK/STAT followed a similar expression pattern, but was Figure 3. Log-centred TMM-normalised FPKM expression levels for genes associated with biogenic amine pathways.Results of Kruskal-Wallis tests are given above each plot facet.Full names and accession numbers of the genes included can be found in Suppl.material 1: Table S3.
only significantly down-regulated in Europe compared to Argentina (p < 0.05; Figure 4a, Suppl.material 1: Table S5).Conversely, genes associated with the Toll pathway were consistently significantly up-regulated in the introduced range (Figure 4e, Suppl.material 1: Table S5).Genes associated with Imd and JNK pathways did not exhibit a clear expression pattern between the native and invaded ranges (Figures 4c-d, Suppl.material 1: Table S5).S4 and full results of statistical tests in Suppl.material 1: Table S5.

Virome analysis
De novo-assembled transcripts that matched virus sequences using BLASTx searches were 102-2091 bases long (mean: 303.5 bases; median: 151; Suppl.material 1: File S4).Virus accumulation curves suggest a sufficient sampling size for most introduced regions, al-though a plateau was not reached in Europe or the native range (Figure 5a).Fewer virus species were detected in California, Australia and New Zealand samples compared to Argentina and Europe and all study populations harboured a core virome of 9 known viruses (Figure 5b, Table 1).These ubiquitous viruses included Linepithema humile picorna-like virus 1 (LhuPiLV1), bunyan-like virus 1 (LhuBLV1), partiti-like virus 1 (Lhu-PLV1), toti-like virus 1 (LhuTLV1), C virus 1 (LhuCV1), rhabdo-like virus 1 (LhuRLV1), polycipivirus 2 (LhuPcV2), Linepithema humile virus 1 (LHUV-1) and Kashmir bee virus (KBV).Virus presence was confirmed with RT-PCR for a subset of viruses (Table 1).Table 1 c Overall viral loads per sample, showing contributions for all detected viruses in the dataset, expressed in TMM-normalised TPM scaled to the Argentine ant library size.Samples are identified with collection site and replicate in subscript (Suppl.material 1: Table S1).
Viral loads, expressed as the sum of viral transcripts detected, were extremely variable at the individual virus species level, as were overall viral loads computed by summing counts for all virus sequences (Figure 5c).We found no significant differences in overall viral loads amongst regions (KW = 3.432, df = 4, p = 0.488).

Discussion
We examined possible functional adaptations underlying the success of a globally invasive pest by investigating transcriptome-wide expression profiles associated with range in the Argentine ant.First, we identified the most differentially expressed genes amongst regions across the introduction pathway, as well as modules of co-expressed genes.We then further investigated gene expression profiles associated with immune pathways, as well as biogenic amine signalling.We also identified viral transcripts present in the libraries to measure viral diversity along the introduction pathway, as well as overall viral loads.Unsupervised multi-dimensional scaling analysis, based on normalised expression of all expressed transcripts, showed range and region-driven clustering.Similarly, hierarchical sample clustering, based on the most differentially expressed genes, showed perfect range-wise and close to perfect region-wise clustering.Functional analysis of differentially expressed genes in both analyses indicated a number of genes associated with a wide range of biological processes, including immunity.Further scrutiny at specific gene groups pointed to consistent range differences in genes associated with biogenic amines and key immune pathways.We also found lower viral Behavioural traits key to the ecological success of ants, such as foraging and aggression, are modulated by biogenic amines (Kamhi and Traniello 2013;Kamhi et al. 2017).We investigated if expression of genes related to the dopaminergic, octopaminergic and serotoninergic neural pathways could be associated with the Argentine ant range, potentially contributing to its invasion success through, for instance, increased foraging activity, increased interspecific aggression and/or lower intraspecific aggression.Interestingly, we did not find behavioural differences in foraging activity between native and invaded ranges, and we found no variation in sensitivity to OA supplementation in the diet provided to the experimental colonies (Felden et al. 2018).Here, genes associated with biogenic amines, including OA, displayed clear range-specific expression pattern, suggesting an association with the introduction process.Biogenic amines are associated with a wide range of behaviour in insects, including key social behaviour in ants (Kamhi et al. 2017) and it is interesting that the vast majority of genes associated with biogenic amines exhibited the same trend.Although the insect brain is comparatively smaller than that of a vertebrate brain, it is a complex organ with largely heterogeneous regions that fulfil different functions and likely exhibit widely variable gene expression patterns.Further analysis of neuropil-specific expression patterns will be useful to elucidate variation in the molecular basis of behaviour.
In our study, all Argentine ant populations harboured a core virome of nine viruses, but we found lower viral diversity in most of the introduced range (i.e.California, Australia and New Zealand) compared to Argentina and Europe, where virus species richness was the highest.Variation in viral diversity has similarly been shown in another widespread invasive ant, the red imported fire ant, Solenopsis invicta (Yang et al. 2010).We found that the RNAi and Imd pathways, as well as the JaK/STAT and JNK pathways, are down-regulated in the introduced range.In honey bees, the JaK/STAT and RNAi immune pathways are primarily directed against viruses, whereas the Toll, Imd and JNK pathways are more ubiquitous and associated with defence against bacteria and other pathogens (Evans et al. 2006;Brutscher et al. 2015Brutscher et al. , 2017;;McMenamin et al. 2018).Interestingly, genes associated with Toll pathway in our results did not show a clear expression pattern (i.e.candidate gene approach), but a number of Toll-associated genes were consistently up-regulated in the introduced range (differential gene expression analysis) and/or present in the gene module (WGC-NA), suggesting a complex pattern of variation in immune responses between native and introduced populations.Changes in bacterial communities in the Argentine ant populations which we studied (Lester et al. 2017) with both acquisition and loss of bacteria may also be related to the variation in Toll genes expression that we observed.Host/pathogen interactions can be unique and much work remains to be done to fully understand specific immune responses associated with different pathogens and symbionts in Argentine ants.In Argentine ants, it has been that different viruses can trigger distinct immune responses, sometimes associated with specific virus families (Lester et al. 2019) and exposure to fewer viruses or a different viral community may result in population-wide differences in immune gene expression.Interactions amongst viruses within individual hosts and host communities appear to be complex (Viljakainen et al. 2018) and it is therefore hard to disentangle the respective weight of viral loads and diversity.Viral diversity may, however, be an informative metric in the context of biological invasions as it is less susceptible to variation at the population level.
Lower viral pressure in the introduced range may promote invasion success in two ways.First, ants are known to harbour a range of viruses (Valles et al. 2007(Valles et al. , 2018;;Sébastien et al. 2015;Cooling et al. 2017;Gruber et al. 2017;Viljakainen et al. 2018;Valles and Rivers 2019) and lower viral pressure may directly increase their fitness.For instance, viral infections in the red imported fire ant have been shown to influence colony foraging performance (Hsu et al. 2018).Second, in line with the EICA framework and supported our results of lower expression of certain immune pathways primarily against viruses, relaxed viral pressure may allow reallocation of resources away from immunity to other functions that increase the invader's fitness.Furthermore, different types of viruses are known to trigger different physiological responses (Lester et al. 2019).Therefore, plastic allocation of resources to specific immune responses with respect to variable viral exposure amongst regions may also increase the ant's competitive abilities and persistence of populations (Lester and Gruber 2016).
A significant proportion of genes up-regulated in the introduced range, compared to that of the native range, were surprisingly associated with vision.Some ant species heavily rely on vision to orientate themselves while foraging (Cheng et al. 2014), but Argentine ants have comparatively smaller eyes and are likely to have under-developed vision.Chemoreception appears to be the primary modality of Argentine ant's sensory biology, as indicated by the large number of chemoreceptor genes in their genome (Smith et al. 2011).Nevertheless, our results may indicate that vision has been overlooked in past studies.Alternatively and partly because gene annotation in non-model organisms is rarely definite, these genes may have been co-opted for other functions in Argentine ants.Throughout our analyses, we also found many genes associated with functions that are difficult to tie to specific biological processes and current resources do not allow for more than speculation on their association with the Argentine ant invasion.Our study does, however, point to several genes and gene families that are worthy of further investigation towards understanding patterns of invasiveness in Argentine ants and other species.
The introduced Argentine ant populations included in this study (i.e.California, Europe, Australia and New Zealand) all belong to the same global 'supercolony' characterised by high genetic similarity and absence of intraspecific aggression, which suggests a common origin (Corin et al. 2007;Van Wilgenburg et al. 2010;Vogel et al. 2010;Suhr et al. 2011).It is possible that the similarity in gene expression patterns, which we observed within the introduced range, is the result of a founder effect at the time of the primary introduction instead of a response to new selective forces associated with the introduced range, as hypothesised in Felden et al. (2018).In order to elucidate the underpinning causes for such variation in gene expression between native and introduced ranges, further investigations should include quantification of genetic variation in functional or regulatory loci related to the signalling pathways discussed herein in these populations, as well as populations from independent introduction events (e.g.Catalonian supercolony in Europe).

Conclusions
We found that native and introduced Argentine ant populations exhibit distinct transcriptomic signatures.Genes associated with biogenic amines were consistently upregulated in the native range, suggesting variation in the molecular basis of behaviour between the native and introduced range.We also observed lower viral diversity in most of the introduced range, which was associated with differential regulation of immune pathways, most notably in the RNAi pathway involved defence against viruses.We provide the first evidence that native and introduced Argentine ant populations are characterised by transcriptomic variation that may reflect region-specific functional adaptations and contribute to the invasion success of the Argentine ant.

Figure 1 .
Figure 1.Gene expression clustering of samples of Argentine ants in their native and introduced ranges, based on multi-dimensional scaling of TMM-normalised counts per million (CPM) with low expressed transcripts filtered out.Euclidian distances between samples are computed from genes with the largest standard deviations.Regions are indicated as Argentina (AR), California (CA), Europe (EU), Australia (AU) and New Zealand (NZ) and collection site shown as subscript (details in Suppl.material 1: TableS1).

Figure 2 .
Figure 2. Details of the set of genes consistently differentially expressed between the Argentine ant native and introduced ranges, showing Z-score scaled log2-transformed TMMs with false discovery rate (FDR) < 0.05 and fold change (|FC|) > 1.1.Trees are based on average Pearson correlation of gene expression.

Figure 4 .
Figure 4. Variation in candidate immune gene expression in the Argentine ant native (Argentina, AR) and introduced (California, CA; Europe, EU; Australia, AU; New Zealand, NZ) ranges.Log-centred TMM-normalised FPKMs for genes associated with the immune pathways (a) RNAi, (b) JaK/STAT, (c) Imd, (d) JNK and (e) Toll.P-values for GLMs or Kruskal-Wallis tests, testing for differences between introduced regions and Argentina are given above each introduced range boxplot, compared to Argentina.Details of the genes, included in the datasets, are found in Suppl.material 1: TableS4and full results of statistical tests in Suppl.material 1: TableS5.

Figure 5 .a
Figure 5. a Virus accumulation curves in Argentina (AR), California (CA), Europe (EU), Australia (AU) and New Zealand (NZ) as detected from RNA-seq data.Error bars indicate standard deviations b Venn diagram showing viral diversity and overlap amongst regions.Detail of the data is given inTable 1c Overall viral loads per sample, showing contributions for all detected viruses in the dataset, expressed in TMM-normalised TPM scaled to the Argentine ant library size.Samples are identified with collection site and replicate in subscript (Suppl.material 1: TableS1).

Table 1 .
Viral diversity detected via RNA-seq, presence is indicated with 1, absence with a dash.Ubiquitous viruses are shown in bold.Virus presence confirmed by RT-PCR is shown with an asterisk.introduced range, highly variable viral loads between samples within and amongst regions, but no difference in overall viral loads amongst regions.