Evidence for functional convergence in genes upregulated by herbivores ingesting plant secondary compounds

Background Nearly 40 years ago, Freeland and Janzen predicted that liver biotransformation enzymes dictated diet selection by herbivores. Despite decades of research on model species and humans, little is known about the biotransformation mechanisms used by mammalian herbivores to metabolize plant secondary compounds (PSCs). We investigated the independent evolution of PSC biotransformation mechanisms by capitalizing on a dramatic diet change event—the dietary inclusion of creosote bush (Larrea tridentata)—that occurred in the recent evolutionary history of two species of woodrats (Neotoma lepida and N. bryanti). Results By comparing gene expression profiles of two populations of woodrats with evolutionary experience to creosote and one population naïve to creosote, we identified genes either induced by a diet containing creosote PSCs or constitutively higher in populations with evolutionary experience of creosote. Although only one detoxification gene (an aldo-keto reductase) was induced by both experienced populations, these populations converged upon functionally equivalent strategies to biotransform the PSCs of creosote bush by constitutively expressing aldehyde and alcohol dehydrogenases, Cytochromes P450s, methyltransferases, glutathione S-transferases and sulfotransferases. The response of the naïve woodrat population to creosote bush was indicative of extreme physiological stress. Conclusions The hepatic detoxification system of mammals is notoriously complex, with hundreds of known biotransformation enzymes. The comparison herein of woodrat taxa that differ in evolutionary and ecological experience with toxins in creosote bush reveals convergence in the overall strategies used by independent species after a historical shift in diet. In addition, remarkably few genes seemed to be important in this dietary shift. The research lays the requisite groundwork for future studies of specific biotransformation pathways used by woodrats to metabolize the toxins in creosote and the evolution of diet switching in woodrats. On a larger level, this work advances our understanding of the mechanisms used by mammalian herbivores to process toxic diets and illustrates the importance of the selective relationship of PSCs in shaping herbivore diversity.


Background
For herbivorous animals, food ingestion can have deleterious consequences [1,2]. At every meal, herbivores confront potential toxins in the form of plant secondary compounds (PSCs). Some herbivores, e.g., the black swallowtail (Papilio polyxenes) and Stephen's woodrat (Neotoma stephensi), have evolved the ability to specialize on plant species with high concentrations of PSCs, whereas other herbivores must consume undefended plants or small quantities of many species of plants with differing PSCs to keep toxin intake below threshold doses [3][4][5][6].
In general, diet selection and diet breadth of herbivores is thought to be governed by the capacity of biotransformation (detoxification) enzymes in the liver [7,8]. Major evolutionary dietary shifts seem to be enabled by an underlying change in biotransformation enzymes [9,10]. Surprisingly few alterations among myriad biotransformation enzymes are required for a change in the capacity to metabolize a new toxin. Insect resistance to pesticides and mammalian resistance to rodenticides can occur through the increased activity of a single biotransformation enzyme [9,11]. Currently, far more is known about the biotransformation mechanisms employed by insect herbivores in a diet shift than their mammalian counterparts [1,9,12].
The vast diversity of biotransformation enzymes presents a challenge in determining which enzymes enable the processing of a particular diet or toxin. For example, more than 300 hepatic biotransformation enzymes have been described for the laboratory rat, Rattus norvegicus (www.reactome.org, [13]). Enzyme activity assays exist for far fewer than the known number of biotransformation enzymes and are not always specific to a gene product [14]. To overcome these challenges, we took a genomic approach and used microarray technology customized for our study species, to identify, on a more global scale, biotransformation gene expression patterns involved in the processing of PSCs. Gene expression is an important phenotypic character, subject to selection from the environment [14][15][16]. Many recent studies have concluded that differences in gene expression among groups are the result of disparate selective forces [17][18][19].
We investigated the independent evolution of biotransformation mechanisms related to diet switching by capitalizing on a dramatic diet change event that occurred in the recent evolutionary history two species of woodrats (Neotoma lepida, N. bryanti). Both species underwent radical dietary changes due to floral shifts during the climate change event of the late Pleistocene and early Holocene (18,700-10,000 years ago). This event resulted in both of these woodrats independently switching from diets of juniper and/or cactus, to one containing creosote bush as it naturally invaded from Mexico and expanded to become the dominant shrub in the major deserts of the U.S. southwest [20,21]. Fossil records of middens of both species indicate that N. bryanti came into contact with creosote 7,000 years before N. lepida [22,23].
The change in diet to creosote from previous plant species represents a marked change in PSC composition. Juniper contains high concentrations of numerous terpenes as well as less abundant tannins [24,25]. The primary PSC in cactus (oxalate) is only degraded by gut microbes, not liver enzymes [26]. In contrast, creosote leaves produce a complex resin containing >300 compounds, primarily aromatic ones, which account for 10-25% of the leaf dry weight [27]. The principal component of resin is nordihydroguaiaretic acid (NDGA), a phenolic lignan with detrimental effects when administered to laboratory rats in doses regularly consumed by woodrats [28][29][30]. The distinct differences in the PSCs of juniper and cactus versus creosote suggest they are metabolized by different biotransformation pathways [14,31]. Comparative data on enzyme activities and gene expression of populations of N. lepida support this contention. Neotoma lepida that feed on juniper have different enzyme activity and gene expression profiles compared to those that feed on creosote even when they are fed the same diet [32][33][34]. Indeed, populations of N. lepida in the Mojave desert have adapted to a diet of creosote as evidenced by their ability to ingest greater quantities of creosote compared to those from the Great Basin desert that feed on juniper and have no prior exposure to creosote [35]. Less is known about the ancestral diet of N. bryanti and no data exist with respect to its ability to metabolize PSCs. However, a recent comparative study on gene expression in N. lepida and N. bryanti on a non-toxic diet revealed similarities in expression profiles in populations of both species that feed on creosote, relative to N. bryanti populations naïve to creosote [36]. These results are consistent with potential convergence in biotransformation strategies of N. lepida and N. bryanti with respect to metabolism of a creosote diet.
To advance our knowledge of the mechanisms involved in diet switching in mammals as well as biotransformation of PSCs in general, we addressed three questions. We asked which genes are induced (upregulated) by woodrats in response to the addition of creosote PSCs to the diet. Gene induction is a conventional pharmacological approach used to identify potentially important enzymes in biotransformation of foreign compounds. Second, we asked whether experience mattered with respect to the biotransformation of creosote. We contrasted the responses of woodrats that historically and currently consume creosote to a population that is ecologically and evolutionarily naïve to creosote compounds. Lastly, we asked whether there was convergence in the pathways used for biotransformation of creosote resin between N. bryanti and N. lepida that regularly fed on creosote bush by comparing the induced transcripts and also by contrasting their gene expression profiles.

Feeding trial
All three populations maintained weight on both the control and creosote diet (Additional file 1: Table S1). However, there was a significant difference across populations with the experienced N. bryanti being heavier than the other two populations.

Microarray quality control
Four of the 24 arrays failed to pass all 9 of Agilent's quality metrics. Three of these passed 8 of 9 metrics, and 1 passed 7 of 9. All were kept in the analysis. On average, less than 0.3% of the features across all arrays were flagged as nonuniform. The clustergram grouping individual woodrats by overall expression profile (all 6286 genes) grouped all the N. lepida in one monophyletic clade and each N. bryanti population in its own monophyletic clade, but did not group individuals by diet within species designations (Additional file 2: Figure S1).

Transcripts induced by creosote
Both species and all populations altered gene expression in response to the creosote diet compared to the control. The experienced N. bryanti had fewer induced (n = 7) and repressed transcripts (n = 7) than either the experienced N. lepida (induced n = 26, repressed n = 20) or the naïve N. bryanti (induced n = 18, repressed n = 9); however, the difference among the groups was not significant (Χ 2 = 1.23 p = 0.54). There were no induced transcripts shared by all groups in this analysis. The experienced and naïve populations of N. bryanti shared two induced transcripts, one of which was for a biotransformation gene (aldo-keto reductase 7A3, Table 1). None of the treatment groups on creosote induced more than five different biotransformation genes ( Table 1). The largest fold induction of a biotransformation transcript was observed in the naïve N. bryanti, which induced sulfotransferase 3A1 by 22×. All other biotransformation transcripts across all groups exhibited inductions of 7-fold or less.

Differential expression between naïve and experienced woodrats
There were significant differences in expression related to experience with creosote between N. bryanti populations. Although the total number of transcripts with greater expression in both the experienced and naïve was similar (199 in experienced vs. 183 naïve, [Table 2A with list of genes in Additional file 3: Table S2]), the experienced population expressed absolutely twice as many transcripts for biotransformation enzymes (Table 2A, Χ 2 = 14.72, df = 1, p < 0.001). Proportionally, 37% of all transcripts expressed to a greater degree by experienced N. bryanti were biotransformation related compared to only 18% of those in the naïve N. bryanti. Many of the differentially regulated biotransformation transcripts coded for the same Genbank accession. The effective gene numbers, calculated from Shannon's H, reflected the same pattern as the total biotransformation transcript counts. For naïve N. bryanti the effective gene number was 19.7 compared to an effective gene number for experienced N. bryanti of 34.8.
The degree to which upregulated genes were expressed was greater, overall, in the experienced N. bryanti (Table 3,  Additional file 3: Table S2). The biotransformation transcript with the greatest expression in naïve N. bryanti and the only one with more than 10-fold higher expression in naïve compared to experienced animals, was a cytochrome P450 (2C65). In contrast, there were three biotransformation transcripts in the experienced N. bryanti that exceeded 10-fold higher expression compared to naïve N. bryanti. These transcripts are related to the functionalization of aldehydes (aldo-keto reductase 1C12) and glucuronic acid and glutathione conjugation pathways (UDP glucuronosyltransferase 2B34, glutathione S-transferase, mu 7).
Significant Gene Ontology (GO) terms (within biological process) were different in the naïve versus experienced N. bryanti. Before comparison, terms across ontology tiers were filtered to include only those with significant, positive z-scores (>2) and ordered by the number of genes in each ontology term. GO terms are reported only if they include 10 or more genes from the array. Naïve N. bryanti had more than twice the number of GO terms that passed these criteria. The top terms related to responses to stimuli; many of the other terms related to stress, damage and cell death (Table 4A). The experienced N. bryanti had fewer significant terms; the function of these related largely to metabolic processes (Table 4B).

Differential expression between experienced woodrats of two species
The experienced N. bryanti had increased expression of more transcripts than N. lepida (134 vs 109) as well as a greater degree to which those genes were expressed (Additional file 4: Table S3). This result was also true with respect to the biotransformation transcripts alone; N. bryanti had greater expression of 3× more transcripts compared to N. lepida on creosote and, in general, these genes were expressed to a much higher degree (Table 2B and Table 5). The effective gene numbers also reflect this pattern. The effective gene number with greater expression in N. bryanti was 20.7, whereas the effective gene number for N. lepida was 14.4. In experienced N. bryanti, transcripts with the greatest expression were related to conjugation of metabolites with glucuronic acid (different glucuronosyltransferases) followed by those related to conjugation with glutathione (glutathione S-transferases; Table 5). In N. lepida, the transcripts with the highest expression were those related to functionalization of aldehydes by aldehyde dehydrogenase (1 L1) and oxidation by Cytochrome P450 (2C65). Compared to N. bryanti, N. lepida also had greater expression of some of the same transcripts that were significantly induced by creosote feeding, i.e., aldehyde dehydrogenase (1 L1) glutathione S-transferase mu-1 (Tables 1 and 5).
Significant GO terms overrepresented by experienced N. bryanti fed creosote were characterized by metabolic and catabolic processes as well as responses to stimuli (Table 6A). GO terms overrepresented by experienced N. lepida are characterized by responses to stimuli and signaling functions (Table 6B).

Discussion
Despite the vast knowledge of drug-metabolizing enzymes in humans and model species, the biotransformation mechanisms used by mammalian herbivores to metabolize PSCs are largely unknown. Until recently, many studies were limited to analysis of one or a few biotransformation enzymes [32,[37][38][39][40]. This study took advantage of microarray technology customized for a unique study system to address this deficit in our understanding. We investigated, on a transcriptomic scale, the genes induced by a particular suite of PSCs (creosote resin) and evaluated whether a common set of genes were expressed by herbivores with varying levels of evolutionary experience with these PSCs. We explored whether experienced herbivores have independently converged on regulation of a similar set of biotransformation genes. We found that biotransformation enzyme expression does indeed vary with ecological and evolutionary experience with creosote and that independent woodrat lineages employ, in part, similar strategies for successfully dealing with these shared PSCs. In addition, the results narrow the field from hundreds of possible biotransformation genes to less than ten candidates. The work provides a testable framework for the changes in expression of biotransformation enzymes that may have occurred as woodrats shifted from one toxic diet to another.

Induced biotransformation genes
We identified a narrow set of candidate genes relevant to the biotransformation of creosote. Surprisingly few biotransformation transcripts were induced by any of the groups fed creosote resin compared to the control diet. Of the hundreds of biotransformation enzymes, only four unique biotransformation transcripts were induced in the naïve N. bryanti, six in the experienced N. lepida, and in the experienced N. bryanti, only a single biotransformation enzyme, an aldo-keto reductase (AKR7A3) was induced ( Table 1). All of these transcripts encode for enzymes that act on substrates similar to the compounds present in creosote, particularly aromatic compounds [27,31]. These enzymes could function in tandem as a pathway to produce the final metabolite excreted in urine and/or feces [41]. For example, epoxide hydrolase acts on aromatic compounds such as naphthalene epoxide whose parent compound, napthalene, is present in creosote [41]. AKR enzymes, in turn, act on the metabolites produced by epoxide hydrolase, whereas glutathione Stransferases and UDP glucuronosyltransferases add conjugates to the metabolites of AKR [42]. Thus, although few  Chi square analyses compare the number of detoxification genes with higher expression to the overall number of genes with higher expression across the two experimental groups.  transcripts were induced in woodrats fed creosote, those induced transcripts produce enzymes that act on substrates similar to those in resin. Furthermore, the candidate genes identified could function in concert with one another in the biotransformation of PSCs in creosote. Few of the transcripts induced by the creosote diet were shared across woodrat groups. The conspecific populations shared induction of AKR7A3. The superfamily of AKR enzymes act on a broad variety of substrates, particularly aldehydes and ketones, and the 7A3 isoform metabolizes aflatoxin B1 [42]. Epoxide hydrolase was induced by the naïve N. bryanti and experienced N. lepida, and also the experienced N. bryanti but in this case, at less than the 2-fold cut-off. This enzyme adds water to epoxides that could otherwise cause toxicity or mutation [31]. Epoxides are often formed during the biotransformation of aromatic hydrocarbons [31], which are common in creosote resin [27]. The shared induction across all three groups of woodrats suggests epoxide hydrolase could be critical to the biotransformation of resin. Lastly, the naïve N. bryanti and experienced N. lepida induced different sulfotransferases. This superfamily has high affinity to myriad substrates and metabolizes compounds (e.g., polyaromatic hydrocarbons) present in creosote resin [31]. Sulfotransferases biotransform xenobiotics by the addition of a sulfur co-factor, generated from dietary cysteine. Because the availability of cysteine is often limited, sulfation capacity is often much lower than other conjugation pathways. It is possible that the woodrat diet contains ample cysteine or its precursors, serine and methionine. Consistent with this idea, both experienced populations upregulated methyltransferases involved in serine and methionine metabolism (Tables 3 and 5). Alternatively, sulfation capacity may be greater in woodrats compared to other mammals. However, previous work suggests the capacity for sulfation is exceeded in both naïve and experienced N. lepida at low levels of creosote ingestion [43]. The role that sulfation plays in the biotransformation of creosote resin deserves further attention.

Alterations in biotransformation mechanisms as a result of diet shifts
The change by woodrats from an ancestral diet of juniper and/or cactus to that of creosote likely required a different set of biotransformation enzymes to process the disparate suites of PSCs in these plants. The naïve and experienced N. bryanti had increased expression of about the same number of genes on creosote compared to the control diet. However, within that pool of genes, the naïve N. bryanti expressed fewer transcripts with biotransformation functions than the experienced animals. A similar pattern has been observed elsewhere; woodrats feeding on a novel diet (i.e., animals naïve to the diet) expressed relatively few biotransformation transcripts and more transcripts related to cellular function compared to feeding on their native diet [33]. The animals on the novel diet were thought to be unable to marshal the appropriate biotransformation response, and the upregulation of genes with diverse cellular functions may have been to prevent or respond to the physiological consequences of the PSCs. Indeed, the functional analyses of the gene ontologies for differentially expressed genes in N. bryanti are consistent with this hypothesis ( Table 4). The transcriptomes of naïve woodrats were enriched in GO classes related to considerable physiological stresses (e.g., death, response to wounding). In contrast, those enriched in experienced N. bryanti were related to the processing of toxins.

Functional convergence
There is evidence for functional convergence between the two species with respect to the transcripts induced on a creosote diet. As discussed above, the transcript for epoxide hydrolase was the only biotransformation transcript significantly induced by both species; however, it was expressed at levels below the 2-fold cut-off used in the initial analysis. The top biotransformation transcript induced by N. lepida (aldehyde dehydrogenase) and the only biotransformation transcript induced by N. bryanti (aldo-keto reductase) both encode for enzymes that metabolize similar substrates, particularly aldehydes. Creosote resin contains at least 300 compounds, including aldehydes and compounds that may have aldehyde functional groups after biotransformation by other enzymes [27,44]. Aldehydes are reactive compounds that cause considerable cellular damage. "List" indicates the number of genes highly expressed within the term; "Gene set" indicates the total number of genes in that ontology included on the array.  Inadequacies in the biotransformation of aldehydes are connected to a number of human diseases [45,46]. Given the potential for damage, aldehydes in creosote are likely a significant selective force resulting in the convergence on the induction of transcripts related to aldehyde metabolism by both species. Further demonstration that the two experienced woodrats exhibit a parallel response to the ingestion of creosote resin is witnessed in their overall gene expression patterns. Many of the biotransformation transcripts that were more highly expressed in the experienced N. bryanti have analogous functions with those that were more highly expressed in N. lepida (Table 5). For example, both species have higher expression of different isoforms of aldehyde dehydrogenases (ALDH8A1 versus ALDH1L1). A similar pattern is present for alcohol dehydrogenases, Cytochromes P450, methyltransferases, glutathione S-transferases, and sulfotransferases. One notable exception to this pattern is that N. lepida did not have a corollary in the UDP-glucuronosyltransferases (UGTs), a superfamily of biotransformation enzymes responsible for glucuronidation. Neotoma bryanti had higher expression of at least 6 different UGTs, one of which was expressed 26 fold higher. The UGTs are considered some of the most versatile of the conjugation enzymes due to their wide range of substrates [31]. Previous studies reported that N. lepida fed creosote increased glucuronidation thereby demonstrating it is indeed an important pathway in N. lepida [32,43]. It is possible that both species use the same UGTs for creosote biotransformation but that N. bryanti utilizes this pathway to a greater extent resulting in the higher gene expression values observed in this study.
The last line of evidence for convergence stems from GO results. The experienced populations fed creosote share three ontologies (response to stimulus, oxidationreduction process, lipid metabolic process), which all relate to biotransformation (Table 6).

Constitutive differences
While many biotransformation enzymes are induced when an animal is exposed to xenobiotics, some are expressed constitutively at high levels [31]. It is difficult to identify such baseline differences in expression using microarrays. For example, a recent study on an herbivorous insect that specializes on glucosinolate-rich plants did not identify any induced biotransformation transcripts when the insects were fed glucosinolates compared to a control diet [47]. Constitutive expression could explain this outcome. In our study, direct comparisons of gene expression between populations fed creosote provide insight into how animals differ in baseline biotransformation, especially when coupled with the results from the induction study. The experienced N. bryanti induce very few transcripts, only one of which has known biotransformation function, relative to both the naïve N. bryanti and the experienced N. lepida. However, direct comparisons of experienced N. bryanti with both the naïve conspecific and experienced congeneric when all are fed creosote, revealed many more differences in expression (Table 2A and B). Overall, the experienced N. bryanti had greater expression of significantly more biotransformation enzymes. This pattern coupled with the paucity of induced transcripts suggests that the experienced N. bryanti may be constitutively expressing biotransformation enzymes at a higher level. It is possible that N. bryanti's longer historic exposure to creosote caused an increase in baseline production of biotransformation enzymes such that they continually express the appropriate combination of enzymes to efficiently detoxify creosote PSCs [22,23]. Moreover, the experienced N. lepida GO terms include many for regulation and signal transduction, suggesting processes in flux, whereas these terms are absent from the experienced N. bryanti results (Table 6).

Conclusions
Herbivores and the plants on which they feed represent one of the most common interactions in nature, yet we know relatively little about the mechanisms that herbivores, particularly mammalian ones, employ to overcome plant toxins [1,48]. The application of recently developed genome-based approaches is enabling researchers to more thoroughly investigate the responses of herbivores to plant defensive compounds. This work capitalized on the availability of a recently developed microarray to identify the genes induced by mammalian herbivores when feeding on PSCs. The results lay the requisite groundwork for future functional and evolutionary studies of the genes involved in the metabolism of creosote toxins and the evolution of diet switching in woodrats.

Trapping and feeding trial
We trapped two species, N. lepida (desert woodrat) and N. bryanti (Bryant's woodrat  "List" indicates the number of genes highly expressed within the term; "Gene set" indicates the total number of genes in that ontology included on the array. (lat: 33.53367, long: −117.54965, Orange Co., CA, USA). This population has no evolutionary or ecological exposure to creosote bush, and is considered "naïve" to creosote PSCs. In preliminary trials, these animals did not maintain body mass on diets of 4% creosote resin. All experimental procedures were approved by the University of Utah's Institutional Animal Care and Use Committee (#07-02015). Eight woodrats from each of the three groups (experienced N. lepida, experienced N. bryanti, naïve N. bryanti) were divided into two treatments of four individuals (n = 2 males, 2 females). One treatment was fed powdered rabbit chow amended with resin from creosote bush; the other group was fed a control diet of the same powdered rabbit chow with no additions. Resin was extracted and creosote diets prepared as in [34].
Animals given the resin diet were fed a gradually increasing amount of resin over an eight-day period to allow induction of biotransformation enzymes [31]. This group received powdered chow with 0% resin for three days, 1% creosote resin (dw/dw) for two days, followed by 2% creosote resin for three days. The final 2% concentration was selected to accommodate the naïve N. bryanti since they were unable to tolerate more resin without excessive mass loss (>10% initial weight). Animals in the control treatments received 0% resin diet for all eight days. Animals had water ad libitum, and diet was presented daily. Animals were weighed daily and body masses were compared within and between groups using a repeated measures ANOVA with time as the repeated variable and population and diet as factors.
Animals were dispatched using CO 2 asphyxiation; liver tissue was preserved in RNAlater (Sigma) and archived at −80°C. RNA was extracted (RNAqueous) for the microarray analysis.

Microarray analysis
Liver samples were analyzed using a custom microarray, built from the hepatic transcriptome of a single desert woodrat [36]. Two primary probe groups were designed from the transcriptome. The target contigs for the first group (n = 943) were woodrat contigs whose annotation matched a list of biotransformation genes extracted from Norway rat arrays that we had previously used in studies with woodrats [34,49]. These probes were all labeled "biotransformation." The target contigs in the second group consisted of contigs from the woodrat transcriptome that had ≥75% sequence identity with a rodent BLAST match in the region of highest homology (n = 5355). All probes were screened for quality and replicated on the array.
The RNA quality and concentration of each sample was evaluated, and the Agilent One-Color Quick Amp Labeling Kit was used to generate fluorescently labeled cRNA for hybridizations. Additional details on microarray processing methods are available in [36]. Microarray features were extracted using Agilent Feature Extraction software version 10.5.1.1. All control spots, non-uniform spots and population outlier spots were removed from the dataset; intensity values were log 2 transformed using Agi-lentFilter, software designed to simplify the processing of Agilent data. Data from duplicate probes were combined, but different probes with the same annotation were maintained separately since it is possible that the original contigs were from different isoforms despite the shared match to a rodent gene. Consequently, many of the resulting lists of differentially expressed genes have multiple seemingly identical entries. For most analyses, we treat these probes as targeting unique genes; in a few specified analyses, we have combined or discarded duplicate probes to evaluate the results as conservatively as possible. The performance of the custom Neotoma expression array was previously validated [36].

Species and diet comparisons
To compare gene expression across diet and species, the data for the Neotoma probes was batch uploaded to Genesifter 3.7. Prior to comparing expression profiles, the consistency of transcriptome response was evaluated by comparing overall gene expression profiles across individuals. After normalizing the intensity data, individuals were clustered by gene expression using all woodrat-derived probes (n = 6286). Clustering parameters were distance: correlation, linkage:average, row centered:by genes.
We created a total of five pairwise comparisons in Genesifter. The specific nature of these comparisons are described in the subsequent paragraphs. In all comparisons, the quality requirement was set to 1, and all individuals from both treatments were required to pass. T-tests were performed with alpha = 0.05, and resulting lists of significantly differentially expressed genes/probes were ordered by fold change with a fold change threshold of greater or equal to 2. This approach is less conservative than implementing a statistical control for the false discovery rate (e.g., Benjamini and Hochberg correction [50]), but fold change ranking is more consistent [51].
To identify induction and repression of genes with respect to ingestion of resin, we compared expression on the creosote diet to the control diet for each of the three groups (experienced N. lepida, experienced N. bryanti, naïve N. bryanti). The lists of significantly differently expressed genes were ordered by direction and fold change. Transcripts that were upregulated on the creosote diet were considered induced by creosote and may be indicative of enzymes important in the biotransformation of creosote. Transcripts that were downregulated were considered repressed. To test for conserved or convergent responses, the lists of induced genes were compared for entries present in more than one woodrat group.
In the fourth comparison, naïve and experienced N. bryanti fed the creosote diet, were compared to identify contrasts in expression with respect to prior experience. Using the gene lists generated from these comparisons, we compared the overall number of transcripts with greater expression in each group as well as the subset of these transcripts characterized by having biotransformation function. The standard gene lists treat probe entries with identical Genbank accession numbers as independent listings. To more conservatively compare the relative response of each group to the creosote diet, we quantified the number of unique Genbank accessions in the biotransformation gene lists and, then calculated the Shannon's H Diversity Index from each list. To calculate Shannon's H, each unique Genbank accession was treated as a "species" and the number of times it appeared in the significantly upregulated list was treated as its "abundance". The values for Shannon's H were then used to calculate an effective number of upregulated genes i.e., the number of "species" expected in comparable communities given their Shannon's H. This index is used to contextualize diversity indices across multiple communities [52]. The effective species value (or effective gene value, in this case) is calculated as exp (H).
To examine function more broadly, we used Gene Ontology (GO). GO is a bioinformatics tool that groups individual genes by the function of their products. GO allows the organization of many individual gene results into fewer functionally-defined categories. There are 3 major classifications: cellular component, molecular function and biological process. Within each of these three, function is categorized more and more specifically. Within our GO results, z-scores were used to determine significantly important associations. Ontology terms with z-scores of >2 are considered to be significantly overrepresented in the results; that is, the genes in that ontology are more likely to be upregulated than expected by chance, given the number features on the array within that particular GO term. Within the biological process classification, two gene ontology lists were generated for transcripts with greater expression in naïve N. bryanti and in experienced N. bryanti. Lastly, experienced animals from both species fed the creosote diet were compared to illuminate contrasts in expression between species with similar experience. This contrast was conducted in the same fashion as the naïve versus experienced N. bryanti.