Prolonged Bat Call Exposure Induces a Broad Transcriptional Response in the Male Fall Armyworm (Spodoptera frugiperda; Lepidoptera: Noctuidae) Brain

Predation risk induces broad behavioral and physiological responses that have traditionally been considered acute and transitory. However, prolonged or frequent exposure to predators and the sensory cues of their presence they broadcast to the environment impact long-term prey physiology and demographics. Though several studies have assessed acute and chronic stress responses in varied taxa, these attempts have often involved a priori expectations of the molecular pathways involved in physiological responses, such as glucocorticoid pathways and neurohormone production in vertebrates. While relatively little is known about physiological and molecular predator-induced stress in insects, many dramatic insect defensive behaviors have evolved to combat selection by predators. For instance, several moth families, such as Noctuidae, include members equipped with tympanic organs that allow the perception of ultrasonic bat calls and facilitate predation avoidance by eliciting evasive aerial flight maneuvers. In this study, we exposed adult male fall armyworm (Spodoptera frugiperda) moths to recorded ultrasonic bat foraging and attack calls for a prolonged period and constructed a de novo transcriptome based on brain tissue from predator cue-exposed relative to control moths kept in silence. Differential expression analysis revealed that 290 transcripts were highly up- or down-regulated among treatment tissues, with many annotating to noteworthy proteins, including a heat shock protein and an antioxidant enzyme involved in cellular stress. Though nearly 50% of differentially expressed transcripts were unannotated, those that were are implied in a broad range of cellular functions within the insect brain, including neurotransmitter metabolism, ionotropic receptor expression, mitochondrial metabolism, heat shock protein activity, antioxidant enzyme activity, actin cytoskeleton dynamics, chromatin binding, methylation, axonal guidance, cilia development, and several signaling pathways. The five most significantly overrepresented Gene Ontology terms included chromatin binding, macromolecular complex binding, glutamate synthase activity, glutamate metabolic process, and glutamate biosynthetic process. As a first assessment of transcriptional responses to ecologically relevant auditory predator cues in the brain of moth prey, this study lays the foundation for examining the influence of these differentially expressed transcripts on insect behavior, physiology, and life history within the framework of predation risk, as observed in ultrasound-sensitive Lepidoptera and other ‘eared’ insects.

Predation risk induces broad behavioral and physiological responses that have traditionally been considered acute and transitory. However, prolonged or frequent exposure to predators and the sensory cues of their presence they broadcast to the environment impact long-term prey physiology and demographics. Though several studies have assessed acute and chronic stress responses in varied taxa, these attempts have often involved a priori expectations of the molecular pathways involved in physiological responses, such as glucocorticoid pathways and neurohormone production in vertebrates. While relatively little is known about physiological and molecular predator-induced stress in insects, many dramatic insect defensive behaviors have evolved to combat selection by predators. For instance, several moth families, such as Noctuidae, include members equipped with tympanic organs that allow the perception of ultrasonic bat calls and facilitate predation avoidance by eliciting evasive aerial flight maneuvers. In this study, we exposed adult male fall armyworm (Spodoptera frugiperda) moths to recorded ultrasonic bat foraging and attack calls for a prolonged period and constructed a de novo transcriptome based on brain tissue from predator cue-exposed relative to control moths kept in silence. Differential expression analysis revealed that 290 transcripts were highly up-or down-regulated among treatment tissues, with many annotating to noteworthy proteins, including a heat shock protein and an antioxidant enzyme involved in cellular stress. Though nearly 50% of differentially expressed transcripts were unannotated, those that were are implied in a broad range of cellular functions within the insect brain, including neurotransmitter metabolism, ionotropic receptor expression, mitochondrial metabolism, heat shock protein activity, antioxidant enzyme activity, actin cytoskeleton dynamics, chromatin binding, methylation, axonal guidance, cilia development, and several signaling pathways. The five most significantly overrepresented Gene Ontology terms included chromatin binding, macromolecular complex binding, glutamate synthase activity, glutamate metabolic process, and glutamate biosynthetic process. As

INTRODUCTION
Predator-induced stress has long fascinated biologists for its integrated, scalable effects on prey physiology, behavior (Slos and Stoks, 2008), and even spatiotemporal population demographics (Clinchy et al., 2013). Though a mechanistic understanding of the physiological responses that are induced by predation related stress in vertebrates has been known for some time, researchers interested in similar responses in invertebrate taxa, such as insects, now seek a similar descriptive model. The study of invertebrate stress responses has a rich history, yet the diversity of molecular components induced by various stressors has thus far stymied most attempts at holistic understanding. Recently, however, Adamo (2010Adamo ( , 2017a demonstrated that the early stages of stress responses in insects are homologous, and likely anciently related, to vertebrate neurotransmitter signaling and downstream neurohormonal activation. The challenge remains, then, in describing the varied taxon-and tissue-specific responses seen in insects and elucidating the mechanisms responsible for inducing them. Often before a predator has even localized its prey, a suite of adaptive behavioral and physiological responses which improve the chances of survival (Endler, 1991) are induced in prey organisms which may be eavesdropping on mechanical, auditory, visual, and chemosensory predation cues (Adamo et al., 2013). For instance, moths and butterflies that are sensitive to ultrasound display startle responses when exposed to synthetic broad frequency ultrasound (Roeder, 1966;Ratcliffe et al., 2008Ratcliffe et al., , 2011ter Hofstede et al., 2011) and recorded bat calls (Acharya and McNeil, 1998;Rydell et al., 2003;Ratcliffe and Fullard, 2005), such as changing the course of flight, ceasing flight, accelerating, performing evasive flight maneuvers (Yack, 2004;Yack et al., 2007;Pfuhl et al., 2015), and/or calling back with jamming ultrasound themselves (Corcoran et al., 2009). Upon exposure to ultrasound, non-flying noctuid moths cease movement while many aerial noctuids exhibit evasive flight maneuvers, such as erratic changes in direction, loops, increases in flight velocity, and even falling to the ground (Surlykke and Miller, 1982). Moreover, when exposed to bat calls, many female and male tympanate moths alter their mating behavior by stopping pheromone release or ceasing flight, respectively (Acharya and McNeil, 1998). These behavioral responses, especially when borne out for an extended period of time, may contribute to patterns of stressor-induced gene regulation in insects that may contribute to reports of moths that display modified fecundity and life history patterns following prolonged exposure to recorded and synthetic bat ultrasound in a laboratory setting (Huang et al., 2003;Zha et al., 2008Zha et al., , 2013. For instance, Plodia interpunctella (Lepidoptera: Pyralidae) exposed to short bursts of ultrasound near their hearing range (approximately 50 kHz) respond by modifying mating behavior (Trematerra and Pavan, 1995) and long-term exposure even affects spermatophore quality and larval numbers by up to 75% (Kirkpatrick and Harein, 1965;Huang et al., 2003) while simultaneously reducing F 1 larval weight and growth rates (Huang et al., 2003;Huang and Subramanyam, 2004). Conversely, long-term exposure to broadband ultrasound in Helicoverpa armigera (Lepidoptera: Noctuidae) significantly increased whole-body acetylcholinesterase activity (Zha et al., 2008), the number of spermatophores per female, and the number of eggs laid (Zha et al., 2013).
In order to maintain internal homeostasis during stressful periods, whether osmotically, metabolically, or otherwise, insects and most other forms of life evolved biomolecular signaling cascades, both intra-and extra-cellularly, that often regulate the expression of stress-related genes (Pauwels et al., 2005;Aruda et al., 2011;Yamaguchi et al., 2012;Roszkowski et al., 2016) and resulting behaviors, including vigilance (Lima, 1990;Kight and Swaddle, 2011) and modified activity patterns (Abramsky et al., 2014). Further, these responses mediate cellular metabolism and the degradative effects of prolonged and persistent stressor exposure, including oxidative damage (Slos and Stoks, 2008;Clinchy et al., 2013), protein misfolding (Fleshner et al., 2004;Even et al., 2012), and organelle turnover (Salvetti et al., 2000;Gesi et al., 2002). However, individual cells can respond to stressful conditions by activating transcriptional pathways that usually produce one or more damage-mitigating antioxidant enzymes or protein folding chaperones, such as the heat shock proteins (Hsps). Though these molecular defenses promote physiological homeostasis in the short-term, prolonged periods of stress clearly influence the life history and fitness of many species. Even though biologists have long recognized the importance of stress hormone signaling for initiating behavioral and physiological defenses to predation, the cellular-and tissuelevel mechanisms by which long-term acclimation to predation risk can influence the life history and fitness of prey species remains unclear, particularly among insects.
In this study, we exposed adult male fall armyworm moths to recorded ultrasonic foraging and attack calls of three insectivorous bat species over an 8-h period to test the influence of an ecologically relevant auditory cue of predation on the cellular physiology of the noctuid brain. The fall armyworm, though a non-model species itself, is in the same family as the corn earworm (Helicoverpa zea), whose annotated reference genome was recently published (Pearce et al., 2017) and whose old world sister species, H. armigera, has long been a prominent subject in insect auditory neuroethology studies for its dramatic neurobehavioral responses to ultrasound. The fall armyworm, and many other tympanic moths, thus make prime candidates for describing the biochemical and cellular responses that have evolved to cope with prolonged predation risk in insects. We hypothesized that a broader transcriptomic response would be induced in the brains of cue-exposed relative to unexposed individuals. Further, we predicted this response might involve transcripts pertaining to the following physiological functions: (1) intracellular secondary messenger systems, (2) antioxidant and Hsp activity, and (3) gene regulation.

MATERIALS AND METHODS
Fall armyworm larvae were purchased from Frontier Agricultural Sciences (Newark, DE, United States) under USDA APHIS PPQ 526 permit (P526P-04080) and were shipped over-night as second and third instar larvae. Upon arrival at the Illinois Natural History Survey, Prairie Research Institute, University of Illinois at Urbana-Champaign in Champaign, IL, United States, larvae were transferred to individual 59 mL (2 oz.) plastic cups filled with 10-15 mL of standard lepidopteran diet and reared in an environmental chamber (Percival Scientific, Perry, IA, United States) at 30 ± 1 • C and 75 ± 5% RH, with a photoperiod of 16 h light/8 h dark. Larvae fed ad libitum on a modified standard larval lepidopteran diet (Sims, 1998;Cohen, 2001;Elvira et al., 2010) prepared every 2 weeks. This diet consisted of 13 g agar, 770 mL distilled water, 31.5 g vitamin-free casein, 24 g sucrose, 27 g wheatgerm, 9 g Wesson's salt mix, 10 g alphacel, 5 mL 4 M potassium hydroxide, 18 g Vanderzant's vitamins, 1.6 g sorbic acid, 1.6 g methyl paraben, 3.2 g ascorbic acid, 0.12 g streptomycin salt, 4 mL wheatgerm oil, and 2 mL 10% formaldehyde. We blended the casein, sucrose, wheatgerm, Wesson's salt mix, alphacel, 220 mL distilled water, and potassium hydroxide on high for 5 min, to which we added 550 mL of mildly boiling distilled, deionized water mixed with agar. We then blended the mixture for another 5 min and allowed it to cool to 60 • C before we added Vanderzant's vitamins, sorbic acid, methyl paraben, ascorbic acid, streptomycin, wheatgerm oil, and formaldehyde and blended for a final 5 min. We poured 10-15 mL of the cooled diet into each 2 oz. rearing cup and allowed them to solidify in a cold-room for at least 30 min.
We then placed a larva into each filled cup and secured a lid in which two holes had been punched using a No. 1 insect pin. Once a larva cleared its gut before pupation, we transferred it to a shallow Tupperware container (29.4 cm × 15.1 cm × 10.5 cm) filled with 3.5 cm of loose potting soil (SunGro Horticulture, Vancouver, BC, Canada). Once per day, this soil was sifted gently by hand to extract any pupae, which were placed in a separate 30.48 cm 3 mesh cage (BioQuip Products, Inc., Compton, CA, United States) with a mesh-size of 51.15 holes/cm 2 within the environmental chamber until emergence.
Upon emergence, adults were transferred to a similar mesh cage and allowed to mate. Twice daily, we saturated the sides of the mesh cage with a 10% sucrose solution to allow feeding. To avoid the possible confounding effects of shipment and the change in diet undergone by the generation of larvae received from Frontier Agricultural Sciences, F 1 eggs were collected daily from within this cage and placed in small plastic containers within the rearing chamber. Once hatched, we reared F 1 larvae as above until emergence as adults.

Predator Cue Exposure
A random sample of four control and four experimental F 1 adult males (sex determined by visual inspection of terminal pupal abdominal segment) were selected for use in trials 24-48 h post-eclosion. Females were not used, as female noctuid moths broadcasting pheromones are often sedentary (Stelinski et al., 2014) and may be preyed upon less frequently by aerial-hawking insectivorous bats. Three individual recordings were sampled at 480 kHz, 16-bit format and concatenated with 10 s of silence between each call. The calls consisted of (1) a 4.27 s Molossus molossus (Chiroptera: Molossidae) attack call, (2) a 1.51 s Myotis nigricans (Chiroptera: Vespertilionidae) foraging call, and (3) a 2.92 Saccopteryx bilineata (Chiroptera: Emballonuridae) foraging call. These three neotropical bat species were selected specifically because the neurophysiological response of S. frugiperda auditory neurons to these species' calls have been explicitly described (Mora et al., 2014), they each represent a ubiquitous species throughout much of S. frugiperda's range in the Americas (Mora et al., 2004;Jung et al., 2007;Surlykke and Kalko, 2008), and they likely represent novel predators for the lab-reared, United Statesbased S. frugiperda colony used in this study. Further, these species produce calls of varying amplitudes and frequencies that together span the known response curve of the S. frugiperda tympanum (Mora et al., 2004(Mora et al., , 2014. Specifically, M. molossus, M. nigricans, and S. bilineata broadcast at 20-50 (Mora et al., 2004), 50-85, and45-55 (Jung et al., 2007) kHz, respectively, whereas S. frugiperda responds optimally to sounds within 20-50 kHz (Mora et al., 2014). The individual sound files were processed in Audacity v. 2.1.0. to reduce background ultrasound by applying a 20-dB noise reduction filter to frequencies lower than 30 kHz with moderate sensitivity (10.0) and re-sampled each file at 195.3125 kHz to meet the limitations of our playback system. This down-sampling attenuated frequencies greater than 75 kHz (Tucker-Davis Technologies, personal communication), but reproduced the bat calls faithfully within the 20-50 kHz optimal hearing range reported for noctuid moths (Fullard, 1988;Norman and Jones, 2000). The resulting 38 s file was then broadcast on a loop for the 8-h duration of each experimental trial while control trials consisted of an identical setup with no sound played whatsoever. Calls were broadcast via a Tucker-Davis Technologies (TDT; Alachua, FL, United States) System 3 amplifier powering an ES1 electrostatic free-field speaker (TDT) that was situated 30 cm from the center of the cage in a soundproof, anechoic chamber at the Beckman Institute, University of Illinois at Urbana-Champaign in Urbana, IL, United States. The RPvdsEx software suite v. 80 (TDT) was used to process and playback the audio file via the TDT RP2.1 processor, ED1 Electrostatic Speaker Driver, and SA1 Stereo Amplifier tandem setup. Each of the four, 8-h replicate exposure and control trials took place on alternating nights in September 2017 from 22:00 to 05:00.

Sample Preparation and Sequencing
Post-exposure, each moth was placed into a 2 mL vial and immediately immersed in liquid nitrogen. After 30 s, the moth was removed from the vial and transferred quickly to a Petri dish on dry ice. After the head was removed, we immersed it in RNAlater stabilization solution (Life Technologies). Upon immersion, scales on the head capsule were removed by scraping with scalpel, and a 1 mm × 1 mm section of cuticle was cut to expose the brain tissue directly to RNAlater. We then dissected the brain from the head capsule, rinsed it with fresh RNAlater solution, placed it in a 2 mL microtube of fresh RNAlater solution, and stored it at 2 • C until all samples had been collected.
RNA was extracted from each brain using a PicoPure RNA Isolation Kit (Arcturus Bioscience). RNA was eluted in 30 µL of RNase-free water and stored at −80 • C until further analysis. Before freezing, 3.5 µL aliquots were removed from each extract and used for RNA quantification via a NanoDrop (Thermo Fisher Scientific) spectrophotometer and a Qubit fluorometer (Life Technologies) using a Qubit RNA HS Assay Kit (Life Technologies). After a 1:10 or 1:15 dilution based on each sample's concentration, we submitted these subsamples to the Functional Genomics Unit of the University of Illinois at Urbana-Champaign's (UIUC) Roy J. Carver Biotechnology Center to confirm RNA quality with a Bioanalyzer RNA 6000 Pico chip (Agilent).
We then submitted each RNA extract to the UIUC Roy J. Carver Biotechnology Center's High-Throughput Sequencing and Genotyping Unit for library preparation and sequencing. Strand-specific cDNA libraries were prepared using an Illumina TruSeq Stranded mRNA Sample Prep Kit (dUTP based) according to manufacturer specifications and quantified by quantitative polymerase chain reaction (qPCR). The eight samples were multiplexed on a single lane of an Illumina 2500 sequencer and the RNA fragments were sequenced using Illumina's HiSeq SBS Sequencing Kit v4 for 101 cycles with a 100 nt paired-end read length.

Raw mRNA Read Preprocessing
Sequence files were demultiplexed with Illumina's bcl2fstq v. 217.1.14 conversion software. To ascertain raw read quality, we used FastQC v. 0.11.2 (Andrews, 2010) with default settings on each set of reads. We then preprocessed the raw reads by performing adapter trimming, quality filtering, and in silico normalization. Adapter trimming and quality filtering was achieved using Trimmomatic v. 0.33 (Bolger et al., 2014) in palindrome mode to search for and remove adapter sequences and low quality bases. To remove redundant reads and improve transcriptome assembly performance, the remaining reads were then digitally normalized to a coverage depth of 50× via the Trinity transcriptome assembly suite v. 2.1.1 (Grabherr et al., 2011;Haas et al., 2013).

De novo Transcriptome Assembly, Annotation, and Quality Assessment
To our knowledge, there is no publicly available annotated reference genome for Spodoptera frugiperda; therefore, we chose to build a de novo transcriptome assembly with the pre-processed reads using the Trinity assembler v. 2.1.1 (Grabherr et al., 2011). We designated the sequence-specific strand orientation to 'reverse-forward' (RF) when possible. The quality of the resulting transcriptome was then assessed using TransRate v. 1.0.1 (Smith-Unna et al., 2016) andBUSCO v. 3 (Simão et al., 2015). We then utilized the Annocript v. 2.0 automated transcriptome annotation algorithm (Musacchia et al., 2015) to complete sequence-similarity searches on each assembled transcript against the National Center for Biotechnology Information (NCBI)'s non-redundant nucleotide database using BLAST+ v. 2.2.30 (Camacho et al., 2009). We selected the UniRef90 protein database (Boutet et al., 2016) to screen for computationally derived protein annotations. Annocript first downloaded the UniRef90 database, stored it in a MySQL v. 7.3 (Oracle Corporation, Redwood City, CA, United States) database, and indexed it for faster searches (Camacho et al., 2009). Annocript carried out BLASTX searches against the UniRef90 database and reported those hits with an e-value < 1e-5. Annocript output a tab-delimited feature map file containing the collated annotation information for each putative assembled transcript.

Read Alignment, Abundance, and Differential Expression Analysis
Following annotation, we indexed the transcriptome in Kallisto (Bray et al., 2016) using the 'kallisto index' command before aligning each sample's reads against the index using the 'kallisto quant' command to select 250 bootstrap replicates each. In R v. 3.5.1 (R Core Team, 2014), we utilized the packages 'edgeR' v. 3.12.1 (Robinson et al., 2009) and'limma' v. 3.26.9 (Ritchie et al., 2015) to import the estimated read counts and perform DE statistical analyses. First, we used the trimmed mean of M-values (TMM) normalization method (Robinson and Oshlack, 2010) to account for small biases in each sample's overall read library size. To filter out transcripts with low or no expression estimates in one or more grouped replicates (Rau et al., 2013), we calculated the counts per million (CPM) mapped reads for each transcript and removed those with a CPM < 1.
We then visually assessed the presence of batch effects in our data by performing principal components analysis (PCA) on log-transformed CPM expression values across each sample using the 'affycoretools' v. 1.42.0 (MacDonald, 2008) package in R. To account for a large amount of expression variation observed between replicate samples (Figure 1A), we used the 'sva' v. 3.18.0 (Leek et al., 2012(Leek et al., , 2010 package to explicitly model three identified surrogate variables as covariates. After adding these covariates to our dataset, we log-transformed all CPM estimates to prepare for linear modeling. We then used the 'limma' package and its 'voom' function (Law et al., 2014) to fit a negative binomial linear model and proceeded to compute pairwise t-statistics, F-statistics, and log-odds of differential expression for each transcript according to exposure type using empirical Bayes (Smyth, 2004). The resulting differentially expressed transcripts were filtered by selecting only those with false discovery rate (FDR)-adjusted p-values < 0.05 and a fold-change > 2 to account for multiple testing bias on p-value significance (Benjamini and Hochberg, 1995;Benjamini and Heller, 2007).
To produce a heatmap of gene expression across the samples, we scaled each transcript's associated fold-change to FIGURE 1 | (A) Principal components plot showing sample clustering based on the first two principal components of variation in log-based counts per million read estimates for both control (C; black circle) and bat-ultrasound exposed (E; gray triangle) Spodoptera frugiperda moths; numbers (1-4) represent replicate samples from each of the control and exposure groups. (B) Principal components plot after surrogate variable analysis was performed to account for unexpected batch effects showing sample clustering based on the first two principal components of variation in log-based counts per million read estimates for both control (C; black circle) and bat-ultrasound exposed (E; gray triangle) adult male Spodoptera frugiperda moths; numbers (1-4) represent replicate samples from each of the control and exposure groups. (C) Transcript expression heatmap detailing the up-(red) and down-(blue) regulation (log 2 FC) of each transcript relative to the mean expression of the control group across bat-ultrasound exposed (E) adult male Spodoptera frugiperda moths; samples (horizontal axis) and transcripts (vertical axis) are clustered according to expression similarity (stacked multicolored bars).
Frontiers in Behavioral Neuroscience | www.frontiersin.org the mean fold-change observed in all transcripts from the control group. To assess similarity in expression between samples, we used a hierarchical clustering method based on a distance matrix compiled by taking the maximal distance between any two expression values in each sample via the 'fastcluster' package v. 1.1.20 (Müllner, 2013) in R. The resultant base dendrogram of similarity between individual transcripts was then used to identify the most appropriate level at which to cluster our transcripts using the R package 'dynamicTreeCut' v. 1.63-1 (Langfelder et al., 2008). We chose to use the 'hybrid' method to first identify large, base clusters following four criteria: (1) each cluster must contain ≥2 transcripts; (2) transcripts that are too distant from a cluster are excluded, even if they occur on the same branch; (3) each preliminary cluster must be distinct from those clusters near to it; and (4) the tips of each preliminary cluster must be tightly connected. Once these clusters were identified, any transcripts not previously assigned were placed in the closest neighboring cluster. Using 'cdbfasta' v. 0.99 1 , we then retrieved the sequences and Gene Ontology (GO) terms associated with these differentially expressed (DE)

Functional Gene Ontology Term Enrichment and KEGG Pathway Analyses
To obtain a broader perspective on the function of our DE transcripts and how they may be related, we tested their associated annotated GO terms for statistically significant over-and under-representation via GO term enrichment analysis. The background set of transcripts we used to test our DE set against included all GO annotations from base transcriptome. Using the 'Biological Networks Gene Ontology' (BiNGO) plugin v. 3.0.3 (Maere et al., 2005) in the Cytoscape platform v. 3.3.0 (Shannon et al., 2003), we tested for both over-and under-representation using a hypergeometric test at an FDR-adjusted p-value < 0.05. Each differentially expressed transcript was also annotated using the automated BlastKOALA (Kanehisa et al., 2016) KEGG pathway webserver and analyzed manually for functional relevance.

Data Availability
The raw sequence reads have been uploaded to the NCBI Sequence Read Archive (SRA) database (accessions: SRR3406020, SRR3406031, SRR3406036, SRR3406052, SRR3406053, SRR3406054, SRR3406055, SRR3406059) and are also available through the BioProject accession PRJNA318819 2 . The transcriptome has been archived to NCBI's Transcriptome Shotgun Assembly database under accession GESP00000000; the version used here is GESP00000000.1. A repository containing R scripts and output files from all analyses downstream of assembly is also hosted on GitHub 3 .

RNA Extraction, Library Preparation, and Read Processing
Each RNA extract was found to produce satisfactory yields, and these were subsequently used in downstream analyses. Total RNA concentration in each sample was generally consistent between NanoDrop and Qubit estimates, the absorbance ratios signified little if any contamination (A260/280 > 2) and the bioanalyzer assay revealed each sample consisted of high-quality RNA with negligible signs of degradation (RIN > 8; Table 1). Our cDNA fragment lengths after library preparation ranged from 80 to 700 bp, with an average of 300 bp. Each sample produced similar numbers of reads, ranging between 28.3 million and 31.1 million. The average quality scores for each base in each sample were ≥33 (phred-33 scaling), allowing us to proceed without sequencing error correction. Preprocessing steps led to less than 0.12% of reads being removed in each sample, and the GC content of the samples ranged from 42 to 45% post-trimming.

De novo Transcriptome Assembly Statistics
Our transcriptome contained a total of 27,734 putative transcript contigs in total, ranging in length from 124 to 38,522 bp with an average of 1,399.9 bp. The contig N50 of our assembly was 2,933 bp and 40.9% of the contigs exceeded 1,000 bp in length. Out of 303 eukaryotic orthologs used as a reference in BUSCO, we identified 290 (95.7%) complete matches, with 234 singlecopy and 56 duplicate hits, along with three fragmented and ten missing orthologs. Annocript annotated 10,367 (37.38%) contigs with reliable protein annotations from significant (evalue < 10 −5 ) BLASTX hits using the UniRef90 database. Of these hits, 97.0% and 93.6% were annotated to species of Insecta and Lepidoptera, respectively. Mapping GO annotations to these hits resulted in 6,476 GO annotations present in the transcriptome, with 4,075 gene products attributed to biological processes, 815 to cellular components, and 1,586 to molecular function. The top GO terms attributed to the largest numbers of transcript contigs included 'integral to membrane' (GO:0016021), 'nucleic acid binding' (GO:0003676), ' ATP binding' (GO:0005524), 'nucleus' (GO:0005634), and 'zinc ion binding' (GO:0008270).

Read Alignment and Abundance Quantification
On average, 36.59% ± 0.573% (95% CI) of reads from each sample mapped to the transcriptome. TMM normalization resulted in normalization factors ranging from 0.915 to 1.104, which we then multiplied by our actual library sizes to find our final effective library sizes. After filtering low and no expression transcripts with <1 CPM, 17,558 out of 27,734 (63.3%) were retained for DE analysis.

Differential Transcript Expression Analysis
Our initial PCA indicated strong, unexpected clustering of samples along the first two principal axes ( Figure 1A), leading us to use surrogate variable analysis in effort to remove potential unaccounted batch effects. We found three significant surrogate variables that we included in our negative binomial regression model as covariates, resulting in clear clustering of samples by experimental group (Figure 1B). Further, we improved our detection of significant DE transcripts at a FDR < 0.05 with ≥2-fold change in expression from 75 to 290 transcripts after including the covariates ( Figure 1C).

Upregulated Genes
Among the top 10 most highly upregulated genes were a X-linked retinitis pigmentosa GTPase regulator (RPGR) homolog and a mitochondrial calcium uniporter protein, though seven genes were unannotated, including the most highly upregulated transcript, with the remaining transcripts annotated by uncharacterized proteins. Genes also had highly variable absolute log 2 -transformed fold changes (log 2 FC) ranging from 1.29 to 11.45. Additional upregulated genes of interest include the regulatory-associated protein of TOR, axin, inositol 1,4 5-triphosphate 5-phosphatase, Hsp 67B2-like isoform X2, glutathione (GSH) S-transferase 2-like, and the rho GTPaseactivating protein.

Gene Ontology Enrichment Analysis and KEGG Pathway Reconstruction
Of the 146 DE annotated transcripts, 102 (69.8%) displayed GO term sequence identity (Figure 2A). GO term enrichment analysis identified 15 overrepresented and 0 underrepresented GO categories in our exposed samples (FDR-adjusted p-value < 0.05; Table 4). Six of these overrepresented GO terms pertained to glutamate metabolism, biosynthesis, and synthase activity, while dicarboxylic acid biosynthesis and metabolism corresponded to two terms, and oxidoreductase, aminoacylase, flavin mononucleotide binding, chromatin binding, and macromolecular complex binding corresponded to one term each. Notably, 14 of these 15 overrepresented GO terms annotated a downregulated transcript while only a single term pertained to an upregulated transcript. Of note is that the majority of transcripts mapping to significantly enriched GO terms occurred as very low or zero transcript count observations in the exposed relative to the control group. All transcripts mapping to chromatin binding-, glutamate-, integrin-, oxidoreductase-, and aminoacylase-related GO terms exhibited this pattern of "all-or-nothing" transcript expression. As the data included considerable noise, the prevalence of this pattern among the differentially expressed GO annotated transcripts may simply be due to these patterns being the only ones strong enough to discern statistically, though their functional relevance in stress physiology requires further investigation. Our BlastKOALA KEGG pathway reconstruction of the 290 DE transcripts recovered 43 (14.8%) with functional annotations, including 37 pertaining to cellular metabolism, six related to genetic information processing, nine that function in cellular signal transduction to environmental stimuli, five related to cell growth and death, two related to glutamatergic and GABAergic 2 | List of differentially upregulated (log 2 -transformed fold change) transcripts recovered from brain tissue mRNA extractions in bat call-exposed Spodoptera frugiperda adult male moths relative to controls, including the most significant (e-value < 1e-5) BLASTX protein annotation from the UniRef90 database and the organism from which the annotation is derived.    synapses, respectively, and one related to neurotrophin signaling in neurons specifically ( Figure 2B).

A Comparison of Predator-Induced Gene Expression Responses in Other Animals
Our results build on a growing body of literature detailing auditory sensory mode and predator-induced shifts in gene expression in vertebrates and invertebrates (Nanda et al., 2008;Leder et al., 2009;Preisser, 2009;Sheriff and Thaler, 2014;Takahashi, 2014;Harris and Carr, 2016;Adamo, 2017a,b). Several studies have focused on describing the gene expression dynamics of large-scale predator-induced morphological changes that occur in organisms displaying predation-related polyphenisms, including multiple species of Daphnia (Schwarzenberger et al., 2009;Spanier et al., 2010;Rozenberg et al., 2015) and the Hokkaido salamander (Hynobius retardatus; Matsunami et al., 2015). Less striking predator-induced changes also have been 3 | List of downregulated (log 2 -transformed fold change) transcripts recovered from brain tissue RNA extractions in bat call-exposed Spodoptera frugiperda adult male moths relative to controls, including the most significant (e-value < 1e-5) BLASTX protein annotation from the UniRef90 database and the organism from which the annotation is derived.       (Sanogo et al., 2011) and an intertidal snail (Chu et al., 2014). Exposure to auditory cues of aerial hawking bats for 8 h resulted in significant transcriptomic responses, as evidenced by the wide-ranging fold-changes (log 2 FC) in transcript expression reported here. In the brains of predator-stressed sticklebacks, low-to-moderate fold-changes ranged from 2 to 6 (log 2 FC; Sanogo et al., 2011), while predator-induced polyphenic Daphnia displayed changes ranging from 2 to 10 (log 2 FC; Rozenberg et al., 2015). Furthermore, the number of DE transcripts found here is comparable to that found in other RNA-seq studies on predator-induced gene expression among invertebrates. For instance, Daphnia pulex exposed to kairomones of predatory phantom midge (Chaoborus) larvae displayed 256 DE transcripts (Rozenberg et al., 2015), while only three transcripts were differentially regulated in the intertidal snail Nucella lapillus when exposed to seawater that flowed first through a chamber holding a predatory crab (Carcinus maenas) feeding on N. lapillus (Chu et al., 2014). Further, the number of DE transcripts from brain tissue after predator exposure can vary strongly based on predator identity, as shown by Matsunami et al. (2015) who found that Hokkaido salamander larvae exposed to predatory dragonfly naiads displayed 605 DE transcripts, while only 103 DE transcripts were found after exposure to predatory tadpoles. One primary difference between past studies of predator-induced transcriptional changes that must be considered when interpreting the results presented here is the time scale at which cues of predation are presented. In the case of predator-induced polyphenisms, exposure length depends highly on organism life history but ranges generally from a few to several days. Though our study assesses the effects of prolonged, frequent exposure to an auditory cue of predation over a single night, it should be noted that this time scale is much shorter than used in most other studies of predator-induced transcription. Clearly, the degree to which prey respond transcriptionally to cues of predation risk can vary broadly across taxa and no clear pattern has yet emerged. However, the ubiquity with which metazoan life responds transcriptionally to these cues of predation begs the detailed description of these gene pathways, their relevance to physiology and life history, and their evolution throughout the tree of life.

Functional Relevance of Differentially Regulated Genes
Furthermore, our results indicate a broad range of functional annotations related to our DE transcripts. For instance, upregulated transcripts coded for proteins related to cellular signaling, Hsp synthesis, antioxidant metabolism, mitochondrial metabolism, oxidoreductase activity, glutamate synthesis, ionotropic receptor activity, gene regulation, ion transport, and cilium assembly. Downregulated transcript annotations also displayed a large degree of functional variability relating to G-coupled protein signaling, cytochrome P450 activity, chromatin-mediated gene regulation, integrin signaling, glutamate biosynthesis, and voltage-dependent ion channels, among others. Several notable transcript upregulations corresponded to unexpected protein annotations, including a mitochondrial calcium uniporter protein (log 2 FC = 9.66), an RPGR homolog (log 2 FC = 9.86), mutant cadherin (log 2 FC = 5.60), mitochondrial choline dehydrogenase (log 2 FC = 6.78), and acyl-coenzyme A synthetase shortchain family member 3 (log 2 FC = 6.22). The mitochondrial calcium uniporter protein acts as a transmembrane transporter for uptake of calcium ions into mitochondria for use during respiration (Marchi and Pinton, 2014) after these ions are mobilized from intracellular stores by inositol triphosphate. Notably, another significantly upregulated gene among exposed individuals was type 2 inositol 1,4,5-triphosphate 5-phosphatase (log 2 FC = 6.06). In humans, this phosphatase hydrolyzes inositol triphosphate and functions as a signal-terminating enzyme, preventing further calcium release (Ross et al., 1991;Contreras et al., 2010).
The second most upregulated transcript codes for a RPGR homolog, a protein usually associated with cilia development in the photoreceptors of vertebrate eyes (Gakovic et al., 2011), although it localizes to other tissues and cell types as well (Khanna et al., 2005). We suggest that RPGR upregulation may be related to increased cilia development and neuronal connections but since its expression has not been studied in insect eyes or other tissues, further conclusions about the function of this protein under predator-stressed conditions in S. frugiperda cannot be made. Because S. frugiperda brains were excised without compromising pigment-storing ommatidial cells, the RPGR expression pattern observed here likely is intrinsic to brain tissue and may be related to neural tissues extending from innervations of the eye. Notably, the entire suite of phototransduction proteins found in the Drosophila visual system is also found to act in the fly's auditory transduction system, with visual rhodopsins serving mechanical transduction and amplification roles in auditory neurons of the Johnston's organ (Pumphrey, 1940).
Another upregulated transcript that may be related to neuronal development encoded a mutant cadherin protein found in humans. Cadherins are calcium-dependent cell-cell adhesion proteins that are integral in nearly every step of neural development in larval Drosophila (Fung et al., 2009), have been implicated in guiding new neuron development contributing to neural plasticity (Edsbagge et al., 2004), and are even involved in hair bundle development in vertebrate ears (Hirano and Takeichi, 2012). As expression of cadherins is usually repressed and localized only to synaptic areas in mature brain tissues (Hirano and Takeichi, 2012), the fact that it is highly upregulated in predator-cue exposed S. frugiperda coupled with RPGR upregulation suggests that neural plasticity and development of new neural connections upon exposure to novel environmental cues may play key roles in functionally responding to auditory predator cues.
Several strongly downregulated transcripts also mapped to unexpected protein annotations, including a 27 kDa hemolymph protein (log 2 FC = −10.40), DNA 6mA demethylase-like isoform X1 (log 2 FC = −7.18), decaprenyl-diphosphate synthase subunit 2 (DDSS2; log 2 FC = −6.80), FH1/FH2 domain-containing protein 3 (FHOD3; log 2 FC = −7.85), and Pb-1 (isoform log 2 FC = −9.10, −8.52, −7.60, −7.39, −6.94, 8.54). The 27 kDa hemolymph protein family consists of proteins found in diverse insect taxa but their function remains unknown. DNA 6mA demethylase is another enzyme correlated with a highly downregulated transcript. Methylation of 6mA has been studied primarily in prokaryotes, where it serves as the primary mechanism for epigenetic signaling via DNA methylation-as opposed to the primary mechanism found in eukaryotes, 5-methylcytosine methylation (Vanyushin et al., 1968). Demethylases associated with 6mA and 5-methylcytosine serve to remove methyl groups from DNA and RNA, affecting the transcription and translation of affected nucleic acid chains. In plants and vertebrates, 6mA methylation both increases and decreases transcription factor binding (Luo et al., 2015), while in Drosophila melanogaster loss of a putative 6mA demethylase resulted in increased transposon expression (Zhang et al., 2013). Notably, a transcript annotated with histone-lysine N-methyltransferase (log 2 FC = −2.54) and five transcript isoforms annotated with Pb-1 were downregulated after predator-cue exposure, although another Pb-1 isoform was also upregulated. These proteins are involved in histone H3 remodeling and binding, respectively (Chandrasekaran and Thompson, 2007;. Although the functional significance of these downregulated genes in the brain of predator-exposed S. frugiperda is unclear, epigenetic mechanisms appear to be induced in some manner. The enzyme DDSS2 catalyzes a reaction to supply decaprenyl diphosphate for use in ubiquinone-10 biosynthesis. Ubiquinone-10 is concentrated in mitochondria, where it acts as a component of the electron transport chain during aerobic cellular respiration (Ernster and Dallner, 1995), although it also is found in many diverse organelles at lower concentrations. In this context, ubiquinone-10 acts as an electron transport enzyme moving electrons from enzyme complexes I and II to III in the electron transport chain, a function only it and vitamin K 2 are able to perform (Bhalerao and Clandinin, 2012). Ubiquinone-10 also serves as an antioxidant due to its weak electron affinity when reduced. In this state, electrons are held so loosely that the molecule readily gives up electrons to oxidized substrates. For instance, within mitochondria, ubiquinone-10 prevents the oxidation of DNA nucleotides during interactions between peroxidase and DNA-bound metal ions (López et al., 2010;Miyamae et al., 2013). Although the down-regulation of DDSS2 does not directly imply that lower levels of ubiquinone-10 were present in predator-cue exposed S. frugiperda, further studies should examine ubiquinone-10 responses to predator exposure. With knowledge of the increased mitochondrial metabolic activity suggested by several upregulated transcripts discussed previously, it is surprising that DDSS2 is downregulated, as a greater need for electron transport substrates and antioxidants with enhanced energy production might be expected. Clearly, there is still much to learn in elucidating the role of DDSS2, and mitochondrial metabolism in general, in the context of predatorinduced stress responses.
Formin homology 1/formin homology 2 domain-containing protein 3 (FHOD3), another protein that mapped to a highly downregulated transcript in the predator-exposed S. frugiperda brain, acts as an actin regulator with a scaffolding function and has been found, in humans, to affect organogenesis, tissue homeostasis, and cancer-cell invasion (Katoh and Katoh, 2004). Actin, a protein that forms microfilaments and constitutes the actin cytoskeleton in all eukaryotic cells, plays a key role in cellular locomotion and shape (Lodish et al., 2000). FHOD family proteins are thought to bind to the growing barbed-end of actin polymers and serve both to deliver new actin monomers and promote actin polymerization, effectively mediating the growth of the actin cytoskeleton (Bechtold et al., 2014). FHOD family proteins are regulated by rho-GTPases, a member of which was downregulated after predator-exposure. Furthermore, actinbinding Lin11, Isl-1, Mec-3 protein 3 and alpha catenin were also down-regulated and act as a scaffold protein (Barrientos et al., 2007) and a cellular linking protein between cadherins and actin-containing filaments (Geoffrey and Robert, 2000;Drees et al., 2005;Yamada et al., 2005), respectively. Considering that a transcript encoding a mutant cadherin was upregulated in predator-exposed brains as well, these patterns suggest that the actin cytoskeleton is affected by predator-exposure and that changes in cellular morphology and motility may be involved.

Overrepresented Gene Ontology Terms and KEGG Pathway Reconstruction in Predator-Stressed Brain Tissue
Relative to the list of DE transcript annotations in this study, the overrepresented GO terms enriched in the brains of S. frugiperda after predator exposure were generally restricted to three biochemical pathways: (1) chromatin and macromolecule binding, (2) glutamate synthesis and metabolism, and (3) aminoacylase activity, although terms related to oxidoreductase activity, flavin mononucleotide binding, and integrin signaling also were overrepresented. To the best of our knowledge, these GO terms have not been implicated in any other study of predator-induced transcription. The small set of GO-annotated DE transcripts identified here limit the statistical detection of subtly over-and under-represented terms; regardless, we found 15 GO terms to be highly significantly overrepresented in our set of annotated DE transcripts relative to the frequency at which these terms were found in our GO-annotated transcriptome (p < 0.0004). Chromatin binding (p < 0.0000) and macromolecular complex binding (p < 0.0000) were the most highly overrepresented GO terms identified both with 7 out of 102 GO-annotated DE transcripts mapped to these terms. The binding of cellular proteins to chromatin can elicit varied cellular responses, such as transcriptional regulation, DNA replication, and chromatin remodeling (Ricke and Bielinsky, 2005). Considering that transcripts mapping to ash1 and Pb-1 protein annotations were also differentially regulated, the presence of these GO terms again implies that epigenetic modifications seem to be induced upon exposure to predator cues.
The set of GO terms pertaining to glutamate synthesis and metabolism included glutamate synthase activity, as well as glutamate biosynthesis and metabolism, glutamine family amino acid biosynthesis, and dicarboxylic acid biosynthesis and metabolism. Glutamate, an amino acid anion derived from its dicarboxylic state, glutamic acid, is used during protein synthesis, but is the most abundant excitatory neurotransmitter in the vertebrate brain (Locatelli, 2005). Although acetylcholine is the primary excitatory neurotransmitter in the insect nervous system (Wnuk et al., 2014), glutamate also plays an excitatory role (Leboulle, 2012) as glutamate immunoreactivity (Sinakevitch et al., 2001) and glutamate-induced ion currents (Cayre et al., 1999) have been observed in insect neurons. Intriguingly, application of glutamate to the mushroom body brain regions of the honeybee, Apis mellifera, facilitates glutamatergic neurotransmission and olfactory learning (Locatelli, 2005), and glutamate-mediated neurotransmission has also been implicated in the visual and tactile (Liang et al., 2012) sensory systems. Notably, one of the strongly upregulated (log 2 FC = 5.70) transcripts we found mapped to a fragment of the ionotropic receptor 75d (IR75d). Benton et al. (2009) found that IR75d and 69 other IR-family proteins carry ionotropic glutamate receptorlike amino acid positions and surmised that IRs may similarly act in chemosensory neuron signaling. Although 21 of these 69 novel IRs showed transcriptional responses to chemical signals in the Drosophila antenna, including IR75d, the remaining 46 displayed no chemosensory-related expression (Benton et al., 2009). Further, the presence of different IR subtypes on a given neuron also influences synaptogenesis, synaptic activity, and experience-dependent neural plasticity in Drosophila (Thomas and Sigrist, 2012). Knowing that biochemical pathways pertaining to glutamate production were altered in the brains of predator-cue exposed S. frugiperda coupled with evidence that IR75d was upregulated post-exposure, we suggest that IR75d and its relatives may be involved in the development and function of auditory mechanosensory neurons.
Manual analysis of the KEGG pathway reconstruction of DE transcripts revealed a variety of interconnected neuronspecific metabolic and signaling cascades that were affected by bat ultrasound exposure, including the mechanistic target of rapamycin (mTOR)/Akt, MAPK, Wnt, prolactin, Hippo, and calcium signaling systems, and associated regulatory responses, such as p53, renin-angiotensin, and NF-κB transcript expression. Notably, recent research on the conserved function of these biochemical pathways in the nervous systems of metazoan taxa across phyla describes the function and biological relevance of these pathways on an organismal scale (Mattson and Camandola, 2001;Lilienbaum and Israe, 2003;Pan, 2007;Lau and Bading, 2009;Tedeschi and Di Giovanni, 2009;Brown et al., 2012;Lin et al., 2012;Graber et al., 2013;Flentke et al., 2014;Mao et al., 2014;Patil et al., 2014;Layden et al., 2016;Guo et al., 2017;Haspula and Clark, 2018). For instance, synaptic glutamate (Sinakevitch et al., 2010;Thomas and Sigrist, 2012;Li et al., 2016), mTOR/Akt (Guo et al., 2017), intracellular calcium (Kaltschmidt et al., 2005;Lau and Bading, 2009), and prolactin signaling (Brown et al., 2012;Belugin et al., 2013), followed by differential p53 and NF-κB transcription (Kaltschmidt et al., 2005;Lau and Bading, 2009) are each implicated in the apoptotic and synaptic-activity mediated induction of neural plasticity, learning, and memory from diverse taxa spanning arthropods to chordates. Clearly, much work remains to divulge how the vast evolutionary divergences inherent between the conserved physiological cellular signaling and gene networks of most, if not all, metazoan taxa correlate with lineage and ecology-specific organismal responses to diverse stressors, including predation risk.

CONCLUSION AND FUTURE DIRECTIONS
This study demonstrates that exposure to ecologically relevant auditory cues of predation risk in S. frugiperda results in varied but strong patterns of up-and down-regulation of a broad range of protein products within the moth brain. The most strongly up-and down-regulated transcripts found in this study correspond to many cellular functions which include mitochondrial metabolism, glutamate synthesis and metabolism, actin cytoskeleton morphology and motion, axon guidance, neural structure, and epigenetic modifications. This is a promising first step in developing a model for the transcriptional impacts of frequent and repeated exposure to bat predation cues in S. frugiperda, which may represent acute and chronic responses of cells to predator-induced stress. Several novel predator-cue induced transcriptional pathways are implicated in these results and present promising opportunities for future research. These broad predator-induced transcriptional responses are characteristic of those found in previous studies, such as in predator-stressed stickleback fish (Sanogo et al., 2011), Daphnia (Rozenberg et al., 2015), and the Hokkaido salamander (Matsunami et al., 2015). Contrary to our expectations, there is little overlap between previously reported responses to predatorinduced stress, such as neuropeptide production and increased antioxidant activity, and the novel predator-induced functional annotations reported here. However, mitoferrin, a solute carrier responsible for iron uptake by red blood cells in vertebrates, was significantly upregulated in the brains of stickleback fish repeatedly exposed to a chemical cue of predation (Sanogo et al., 2011), although it was downregulated (log 2 FC = −6.06) in S. frugiperda post-exposure. In insects, the function of mitoferrin is less well understood, though D. melanogaster with mitoferrin mutations experienced problems with spermatogenesis and development to adulthood (Metzendorf and Lind, 2010). Apart from this similarity, the novel transcriptional responses to predation in S. frugiperda observed here may be specialized to auditory perception or found only in Lepidoptera. Furthermore, although efforts were made to avoid auditory habituation in this study, expression profiles described here bear similarities to past studies of bird-song habituation in the brains of zebra finches (Taeniopygia guttata), with both resulting in the downregulation of genes pertaining to cytoskeletal dynamics and mitochondrial metabolism (Dong et al., 2009).
One primary limitation of our study is a lack of timeseries expression data that would have bolstered our ability to infer the functional relevance of specific transcripts for both short-and long-term physiological acclimations to auditory cues of predation. Further work, such as comparing expression profiles through time and between frequent and infrequent cue exposures, would aid in parsing the effects due to neural habituation/auditory stimulation, per se, and those related specifically to predator cue exposure. Specifically, producing a detailed time-course transcriptional profile of tissue-specific prey physiology beginning after the first moments of predator-cue exposure and proceeding over the course of hours to days in cue-exposed S. frugiperda or other predator-prey systems would provide comparative insights into the temporal dynamics of stress-induced transcription during acute relative to prolonged exposure to predation risk. Another limitation of this study is a lack of transcript validation via quantitative reversetranscriptase polymerase chain reaction assays, yet we argue the novelty of the system and the foundational datasets we have produced that can inform future hypotheses warrant their use by the scientific community. Another confounding factor that may have contributed to the relatively noisy patterns of expression in exposed S. frugiperda brains shown here is the type of auditory stimulus we used. For instance, by using three bat calls from three different species, we have endeavored to replicate an ecologically relevant cue of predation risk, yet the nightly soundscape a moth is exposed to in situ varies hourto-hour and night-to-night in sound intensity, conspecific and interspecific composition, and many other attributes that we did not incorporate into our experiments. We encourage future investigators to develop high quality, ultrasonic soundscape recordings in relevant field settings ahead of time, when possible, and replicate these via nightly broadcasts of each night's recording. In conclusion, as more diverse, annotated insect genomes become available and the function of more genes are elucidated by experimental and comparative evidence, studies that assess the physiological effects of prolonged predation risk on prey across the tree of life will continue to divulge remarkably conserved patterns of stress-induced molecular mechanisms between lineages.