|Home | About | Journals | Submit | Contact Us | Français|
Toxicogenomic studies are increasingly used to uncover potential biomarkers of adverse health events, enrich chemical risk assessment, and to facilitate proper identification and treatment of persons susceptible to toxicity. Current approaches to biomarker discovery through gene expression profiling usually utilize a single or few strains of rodents, limiting the ability to detect biomarkers that may represent the wide range of toxicity responses typically observed in genetically heterogeneous human populations. To enhance the utility of animal models to detect response biomarkers for genetically diverse populations, we used a laboratory mouse strain diversity panel. Specifically, mice from 36 inbred strains derived from Mus mus musculus, Mus mus castaneous, and Mus mus domesticus origins were treated with a model hepatotoxic agent, acetaminophen (300 mg/kg, ig). Gene expression profiling was performed on liver tissue collected at 24 h after dosing. We identified 26 population-wide biomarkers of response to acetaminophen hepatotoxicity in which the changes in gene expression were significant across treatment and liver necrosis score but not significant for individual mouse strains. Importantly, most of these biomarker genes are part of the intracellular signaling involved in hepatocyte death and include genes previously associated with acetaminophen-induced hepatotoxicity, such as cyclin-dependent kinase inhibitor 1A (p21) and interleukin 6 signal transducer (Il6st), and genes not previously associated with acetaminophen, such as oncostatin M receptor (Osmr) and MLX interacting protein like (Mlxipl). Our data demonstrate that a multistrain approach may provide utility for understanding genotype-independent toxicity responses and facilitate identification of novel targets of therapeutic intervention.
Biological monitoring to assess potential toxicity of chemical and pharmaceutical compounds relies heavily on the availability of sensitive, specific, and widely applicable biomarkers of toxic effects (International Programme on Chemical Safety, 1993). Toxicogenomics has been used at all stages of chemical risk assessment, and it is thought that gene expression changes may be utilized as biomarkers of adverse effects (Casciano and Woodcock, 2006). Current approaches often attempt to classify compounds with the goals of predicting adverse responses to specific chemical classes (Fostel, 2007), understanding the underlying biological mechanism of toxicity (Dix et al., 2006), or identifying key nodes in the toxicity pathway that may serve as effect biomarkers (Fry et al., 2007). Extensive proprietary (Castle et al., 2002; Ganter et al., 2008) and public (Mattingly et al., 2006; Waters et al., 2008) databases containing gene expression profiles and pathological end points derived from rodent and human tissues exposed to a variety of chemicals have been developed, thereby allowing the scientific community to mine the data for toxicity biomarkers of interest.
Many biomarkers of toxicity may be surrogate measures for the genetics of an individual, which can play a major role in determining the threshold of toxicity of a given compound (Lanfear and McLeod, 2007). Compelling research has led to the identification of gene variants that correlate with drug toxicity (Feero et al., 2008), and recent pharmacogenomic research efforts have made significant advances in connecting variability in responses to drug efficacy and/or toxicity to genetic polymorphisms (Weiss et al., 2008). While major research efforts are seeking genetic and genomic markers that could identify individuals susceptible to toxicity, less attention is given to the fact that interindividual variability in responses and genetic control of gene expression may present a challenge for finding robust population-wide expression biomarkers of toxicity responses (Gatti et al., 2007). Indeed, while toxicogenomics has been used widely for the study of toxicity biomarkers across compounds and across species, its usefulness in determining biomarkers that are relatable to a diagnostic toxicity within a genetically diverse human population is limited by a lack of intraspecies comparisons.
To address the need for a biomarker identification strategy that is independent of population heterogeneity, we utilized a mouse Laboratory Strain Diversity Panel (Bogue and Grubb, 2004). The use of a genetically defined panel of mice has advantages over classical toxicology testing strategies that utilize a single inbred or outbred strain because it takes advantage of the vast genetic diversity that is available among inbred mouse lines (Roberts et al., 2007). We hypothesized that toxicity responses across a panel of strains will produce a range of effects expected in a human population and that this phenotypic diversity can be used to identify population-dependent and -independent mRNA transcript biomarkers of response. To test this hypothesis, we selected the model hepatotoxic agent, acetaminophen. We observed a dramatic gradient of acute hepatotoxicity across strains, and the analysis of liver gene expression data revealed 26 genes that correlated with liver necrosis outcome and were not affected by genetic differences between individual strains. Thus, these genes, the majority of which are tightly linked in a cell death and proliferation network, can serve as response biomarkers for acetaminophen-induced toxicity responses across a genetically heterogeneous population.
Male mice (aged 7–9 weeks) were obtained from the Jackson Laboratory (Bar Harbor, ME) and housed in polycarbonate cages on Sani-Chips irradiated hardwood bedding (P. J. Murphy Forest Products Corp., Montville, NJ). Animals were fed NTP-2000 wafer diet (Zeigler Brothers, Inc., Gardners, PA) and water ad libitum and maintained on a 12-h light-dark cycle. Mice utilized in this study comprise 36 inbred strains that are priority strains for the Mouse Phenome Project (Bogue and Grubb, 2004): 129S1/SvImJ, A/J, AKR/J, BALB/cByJ, BTBR T+ tf/J, BUB/BnJ, C3H/HeJ, C57BL/10J, C57BL/6J, C57BLKS/J, C57BR/CdJ, C57L/J, CAST/EiJ, CBA/J, CZECHII/EiJ, DBA/2J, FVB/NJ, JF1/Ms, KK/HlJ, LP/J, MA/MyJ, MSM/Ms, NOD/ShiLtJ (formerly NOD/LtJ), NON/LtJ, NZO/H1LtJ, NZW/LacJ, P/J, PERA/EiJ, PL/J, PWD/PhJ, RIIIS/J, SEA/GnJ, SJL/J, SM/J, SWR/J, and WSB/EiJ. F1 hybrid mice, B6C3F1/J, were also used for phenotypic measurements. These studies were conducted under a protocol approved by the Institutional Animal Care and Use Committee at the University of North Carolina at Chapel Hill.
Mice were singly housed and fasted 18 h prior to intragastric dosing with acetaminophen (99% pure; Sigma-Aldrich, St Louis, MO; N = 3–4 per strain) or vehicle (0.5% methyl 2-hydroxyethyl cellulose; Sigma-Aldrich; N = 2 per strain, except for strains PERA/EiJ, SWR/J, and CZECHII/EiJ (N = 3), as well as strains AKR/J and CAST/EiJ (N = 1, i.e., sufficient tissue was not available). The dose of 300 mg/kg was delivered in 10 ml/kg of vehicle. Dosing was performed at the same time of day (9:00 A.M.) throughout the study as diurnal effects have been shown to affect gene expression in rodent studies (Boorman et al., 2005). Feed was returned 3 h after dosing; animals were necropsied 24 h after treatment (nembutal 100 mg/kg, ip; Abbott Laboratories, Chicago, IL). Livers were quickly excised following exsanguination, and sections of the left lateral lobe were placed in 10% phosphate-buffered formalin for immunohistochemical analyses. Tissues were stored in formalin solution for 72 h prior to transferring tissue to 70% ethanol. Formalin-fixed liver tissue was then embedded in paraffin. Remaining liver from the left lobe was snap frozen in liquid nitrogen and stored at − 80°C for RNA extraction.
Paraffin-embedded liver tissue was cut to 5-μm sections in duplicate and stained with hematoxylin and eosin. Liver injury in the left liver lobe was blindly scored by A.H.H. and confirmed by a certified veterinary pathologist. Necrosis was quantified by unbiased stereology using a point-counting technique (Mouton, 2002). Briefly, a grid with 100 evenly spaced points was overlaid on printed images of liver sections taken at ×100 magnification. The total number of points lying in an area of necrosis was divided by the total number of points lying completely within the entire tissue section to determine a percent necrosis score (0–100%). The necrosis score for each animal in the study is publicly available from the Mouse Phenome Database (http://phenome.jax.org/pub-cgi/phenome/mpdcgi?rtn = projects/details&sym = Threadgill1).
To eliminate variability in transcript expression that might arise between liver lobes, the left liver lobe was selected for the remainder of the data analysis and gene expression profiling. RNA was extracted from the 30 mg of tissue derived from the left lobe of sample livers using the Qiagen RNeasy kit (Qiagen, Valencia, CA). RNA concentrations were measured using a NanoDrop ND-1000 spectrophotometer (NanoDrop Technologies, Wilmington, DE), and quality was verified using the Agilent Bio-Analyzer (Agilent Technologies, Santa Clara, CA). RNA was determined to be of good quality for use in microarray hybridizations if the 28S:16S rRNA ratio was greater than 1.8 and the 260/280 nm absorbance ratio was in the range of 1.9–2.1.
In this study, all RNA samples were hybridized to arrays individually; none were pooled. RNA amplifications and labeling were performed using Low RNA Input Linear Amplification kits (Agilent Technologies). For hybridization, 750 ng of total RNA from each mouse liver was amplified and labeled with fluorescent dye (Cy5). In parallel, 750 ng of a common reference RNA (Icoria, Inc., RTP, NC) was labeled with the fluorescent dye, Cy3, in order to standardize analysis of global gene expression between mouse strains (Bammler et al., 2005). Labeled cRNA was then processed and hybridized to Agilent Mouse Toxicology Arrays (catalog# 4121A, 22,575 features) according to the manufacturer's protocol. Details regarding the microarray probe set on the 4121A array are available via the Gene Expression Omnibus (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc = GPL891). Following hybridization, arrays were washed using a custom protocol developed by Icoria, Inc. Briefly, array gaskets were removed under immersion in wash solution 1 (6× sodium chloride/sodium phosphate/EDTA [SSPE], 0.005% N-lauroylsarcosine). Arrays were washed with wash solution 1 and incubated for 1 min with gentle agitation on a magnetic stir plate. A second incubation was performed in wash solution 2 (0.06× SSPE, 0.005% N-lauroylsarcosine).
Raw microarray intensity values were obtained from Agilent Feature Extraction software (v8.5) and archived in the UNC Microarray Database (http://genome.unc.edu). Raw data are available to the public through this database. The log2 ratio of Cy5/Cy3 intensity was normalized using locally weighted scatterplot smoothing to eliminate intensity bias of features. Transcripts with fewer than 70% available data across samples were excluded from the analysis, reducing the probe list to 15,509 transcript probes. Available data are defined as those probes that are neither saturated nor below the limit of quantification. Intensity ratios were transformed to eliminate hybridization batch effects using the batch normalization feature in Partek Genomics Suite (Partek, Inc., St Louis, MO). Analysis of significant transcripts was performed using an analysis of covariance (ANCOVA) model in Partek in which the main effects were mouse strain, treatment, the interaction of mouse strain and treatment, and the sample necrosis score. Transcripts were called significantly different if the p value was less than a threshold determined by a step-down false discovery rate (FDR, Benjamini and Liu, 1999) (α = 0.01) to correct for multiple comparisons across array features. Heat maps were generated using hierarchical agglomerative clustering.
There is the potential that strain-specific gene expression differences found using microarrays could be produced by single-nucleotide polymorphisms (SNPs) that occur within probe sequences (Alberts et al., 2007). The genomic locations of the probes on the Agilent G4121A microarray were obtained from Agilent Technologies. High-density mouse SNP data containing 7.87 × 106 SNPs were obtained from Szatkiewicz et al. (2008) for the 36 inbred strains used in this study. We found that 4,091 of the 20,868 probes on the Agilent 4121A array contained at least one SNP. For each probe that contained an SNP, we performed Student's t-tests between the C57BL/6J allele-containing strains and those with the opposite allele. Of these, 948 probes contained an SNP for a single strain, providing no meaningful way to carry out a t-test. We performed Fisher's exact test on the remaining 3,143 probes to determine if there were more probes with C57BL/6J allele high expression on the list of significant probes than would be expected by chance. We found a significant effect of probe SNPs within the data set (P = 0.01) and identified 49 probes that were significantly affected at a 5% FDR. We removed these probes from further analysis.
Ingenuity Pathway Analysis (IPA, IPA v. 7.1; Ingenuity Systems, Redwood City, CA) was used to determine canonical pathways that are enriched by the significant transcripts identified by the ANCOVA model for each factor. Significance values were calculated based upon a right-tailed Fisher's exact test that determines whether a pathway is overrepresented by calculating whether the genes in a given pathway are enriched within the data set compared to all genes on the array in the same pathway; p < 0.05 was selected as the cutoff for significance based on IPA threshold recommendations. Only those pathways with a p value above the threshold and having more than two representative genes in the data set were considered significant. The gene network of the 26 response biomarkers was prepared by determining connecting nodes, interactions, and cellular compartments using the IPA software.
At 24 h after dosing with 300 mg/kg of acetaminophen (ig), we observed centrilobular necrosis in the liver consistent with that previously reported for acute doses of acetaminophen (Hinson et al., 2004; James et al., 2003). Necrosis was accompanied by minor inflammatory infiltration into the hepatic parenchyma, and in varying degrees, hemorrhage was also present. Quantitative liver necrosis scores reflective of the proportion of the affected area were obtained from the left liver lobe (Foley et al., 2006) and demonstrated a wide range of toxicity across the panel of inbred mouse lines (Fig. 1). The rank order of sensitivity to acetaminophen-induced liver injury across strains shows that the majority of tested strains (30/36) sustained less than 40% liver necrosis, while 6 strains sustained liver necrosis of between 40 and 100%.
Gene expression values were collected on individual animals in this study (vehicle and acetaminophen-treated mice) and used for principal components analysis to visually examine the patterns in global mRNA transcript differences (Fig. 2). The unsupervised analysis displayed separation of the samples by both treatment and by the amount of liver necrosis sustained in the animal, indicating that gene signatures may be determined that are correlative with liver toxicity due to acetaminophen.
To determine those transcripts in which expression was significantly differentiated among the experimental factors, an ANCOVA model was used. Covariate factors for each individual mouse included the strain (genotype), treatment (vehicle or acetaminophen), the interaction between strain and treatment because of anticipated genotype-specific effects on acetaminophen metabolism and transport, and the liver necrosis score. The number of transcripts significantly changed among each experimental factor is depicted in a Venn diagram (Fig. 3A). Interestingly, the majority of genes (1,511) found to be significantly different between samples in the ANCOVA analysis were attributed to the strain effect, not acetaminophen treatment, or the degree of liver necrosis. These genes were found not to have a significant effect of probe SNPs within the data set (see “Materials and Methods” section). This strain-specific gene set best represents those genes that differ in basal levels among the panel of inbred mouse strains and whose expression is likely to be influenced by genetic polymorphisms (Gatti et al., 2007). Similarly, those genes that were significant for all three factors (strain, treatment, and necrosis) best represent the genes that could yield important information on the mechanism of acetaminophen toxicity but would make a poor biomarker because basal levels are affected by individual genotype.
Next, functional pathway analysis was performed in order to determine biological pathways most affected by the main experimental factors of strain, treatment, or liver injury (Table 1). Several pathways were identified as significant within the necrosis-specific gene set. As expected, those pathways that were most significantly enriched are associated with cellular growth processes, including signaling mediated by interleukin (IL)-10, IL-6, extracellular signal-regulated kinase (ERK)/mitogen activated protein kinase (MAPK), and nuclear factor kappa (NF- κB).
Interestingly, a diversity of canonical pathways was enriched within the gene set that was significant for a main effect of strain alone. These pathways include genes involved in mitochondrial dysfunction as well as metabolic pathways. Interestingly, the pathway for lipopolysaccharide (LPS)/IL-1–mediated inhibition of retinoid X receptor (RXR) function was significantly enriched and included critical mediators such as peroxisome proliferator–activated receptor alpha and tumor necrosis factor. Taken together, the large number of genes and canonical pathways that are affected by strain-dependent gene expression indicate that choice of rodent strain in toxicity risk assessment is of importance.
There were 26 transcripts whose expression was affected significantly by both treatment and by the toxicity outcome (i.e., liver necrosis), but not the subject's genotype (Table 2). We reason that these genes could serve as population-based biomarkers of response. To visualize gene expression changes of the biomarker transcripts across individuals, a heat map was generated (Fig. 3B). A clear gradient of expression changes can be observed for each of these genes depending on the amount of necrosis sustained by an individual mouse. Expression of 17 of these transcripts increased, while 9 genes decreased as liver necrosis increased in acetaminophen-treated mice.
In order to determine whether molecular interactions exist among the population-based transcript biomarkers, a pathway map was constructed using IPA. This analysis revealed that 16 of the 26 population-based response biomarkers are closely linked in a cell death and proliferation network centered on cell cycle regulating genes Trp53, Myc, Jun, and Cdkn1a (p21) (Fig. 4). Closely associated with this network were the cytokine-responsive genes interleukin 6 signal transducer (Il6st) and oncostatin M receptor (Osmr) (Figs. 3C and 3D), as well as the glucose-responsive transcription factor MLX interacting protein like (Mlxipl) and CDC14 cell division cycle 14 homolog B (Cdc14b) (Figs. 3E and 3F).
Decades of mechanistic investigations into the liver toxicity of acetaminophen have concluded that (1) metabolic activation to the reactive metabolite N-acetyl-p-benzoquinone imine and its binding to cellular proteins is an essential initiating event for the toxicity, (2) intracellular events involved in cell death such as mitochondrial dysfunction and formation of reactive oxygen and nitrogen species propagate the injury, and (iii) inflammatory response to cell death in the liver may exacerbate the damage (Jaeschke and Bajt, 2006; Kaplowitz, 2005). Thus, the fact that our study not only identified 26 biomarker genes in which expression across strains was associated with the level of liver necrosis but also showed that 16 of the 26 response biomarker genes are involved in cell death pathways and form a closely linked molecular network confirms a central role for intracellular cell signaling in acetaminophen-induced liver toxicity.
Not only are cell death–related genes mechanistic biomarkers of response across genetically diverse individuals as identified in our work, they also have been shown to be consistently affected and significantly correlated with the acetaminophen-induced liver toxicity phenotype in a multicenter toxicogenomic study (Beyer et al., 2007). The study, conducted at seven different laboratories around the United States, used only one inbred strain, C57BL/6J; however, it showed that Myc is induced by acetaminophen and that a myelocytomatosis oncogene-centered cell death pathway is the most significant network of proteins associated with liver injury in the mouse at 6, 12, and 24 h after treatment with a dose identical to that used in our work. Furthermore, expression of cyclin-dependent kinase inhibitor p21 (Cdkn1a), a central gene in the biomarker gene network, has been shown previously to be required for liver necrosis in rodents (Kwon et al., 2003). In addition, decreased levels of Cdc14b are consistent with increased activation of Trp53 (Kwon et al., 2003), which may be a compensatory mechanism to signal for an increase in cellular repair following acetaminophen overdose. Collectively, 16 of the 26 response biomarker genes identified in our study may be mechanism-relevant biomarkers of liver necrosis that have potential to be used to profile toxicity across individuals and in multiple independent microarray studies.
Importantly, the genes identified in this study are interesting not only as potential biomarkers but also as mediators of acetaminophen-induced cell death and regeneration in liver. For example, the role of OSMR in acetaminophen-induced liver injury deserves attention because genes coding for its two subunits, Osmr and Il6st, were both identified as genotype-independent biomarkers of the liver toxicity outcome. It is known that Il6st expression is essential for the control of the hepatic acute-phase response during liver regeneration (Streetz et al., 2003; Wuestefeld et al., 2003). However, while IL-6 represents one of the best-studied cytokines, there is relatively little known about the biological activities of oncostatin M (OSM), a cytokine secreted by activated T lymphocytes, macrophages, and neutrophils. OSM may have a profibrotic role in liver injury owing to its ability to induce tissue inhibitor of metalloproteinase (TIMP) 1 (Richards et al., 1993) and TIMP3 (Li and Zafarullah, 1998). While OSM has been shown to be increased following acetaminophen-induced liver injury (Masubuchi et al., 2003), Osmr transcript levels have not been shown previously to correlate with liver necrosis end points. Additionally, knockout mice deficient for Osmr display defects in liver regeneration following carbon tetrachloride exposure (Nakamura et al., 2004); more importantly, administration of exogenous OSM ameliorated liver injury in wild-type mice (Nakamura et al., 2004).
In addition, expression of Mlxipl, also known as carbohydrate response element–binding protein (Chrebp), a transcription factor that plays a central role in the dietary regulation of hepatic gene expression by glucose, was decreased as the degree of liver necrosis increased in animals treated with acetaminophen. Several recent studies demonstrated that acetaminophen can affect blood glucose levels (Kendig et al., 2008) and improve glucose tolerance in mice fed a high-fat diet (Shertzer et al., 2008). The former study showed that daily administration of acetaminophen prevented approximately 70% of weight gain compared to mice fed the high-fat diet alone, even at a daily dose that was lower than half of the maximum recommended weight-adjusted human dose (Kendig et al., 2008). In addition, decreases in liver glucose and increases in lipid content were observed in the mouse liver after acetaminophen overdose using nuclear magnetic resonance-based metabolomics (Coen et al., 2004) and may explain the dramatic decrease in Mlxipl transcript levels observed in our work. While further studies need to be conducted to link effects on glucose modulation at subacute doses of acetaminophen with the acute toxic doses used in our study, changes in Mlxipl expression may yield insight into the mechanism of these phenomena.
An important limitation of the animal studies of toxicity mechanisms is the ability to translate the data to clinical findings. A recent study that compared acetaminophen toxicity in the rat and humans showed that, in the rat, gene expression data from peripheral blood cells can provide valuable information about overtly toxic exposure levels (Bushel et al., 2007). Furthermore, based on the subset of 66 genes that were retrieved from a rat blood training set, it was possible to distinguish humans who overdosed on acetaminophen (five cases) from normal individuals (three controls). None of the 26 genes identified in our multistrain study could be matched to human blood transcriptome and serum alanine aminotransferase data from subjects overdosing on acetaminophen reported by Bushel et al. (2007). However, it should be noted that the small sample size of the human data set, as well as the unavoidable variability in the timing of the collection of human blood samples (2 or 5 days after ingestion of acetaminophen), could have been the major factors for lack of mouse-to-human overlap.
The use of toxicogenomics as a tool in toxicology calls for the careful evaluation of study designs. Because one of the major applications of toxicogenomics is to discover biomarkers of toxicity that are relevant to humans, great care must be taken in choosing the appropriate model systems. Traditional risk assessment practices using animal models allow for the control of many experimental factors except for genetics. Although rodent models have been widely used for toxicity testing, their utility is often limited by (1) inaccurate generalizations from a single genome, (2) inability to distinguish small and biologically important changes from background variation, (3) ineffective exploitation of reproducible genetic variation to dissect differential response to chemical exposure, and (4) inefficient use of defined genetic backgrounds to model particular phenotypic profiles observed in human populations.
To address these important limitations, panels of genetically defined organisms, such as inbred mouse lines, that provide a fixed genotype within a particular strain but encompass great genetic diversity across strains are being used more frequently in biomedical research (Festing, 2001). Genetic variation among individuals is reflected in variations in gene expression levels (Schadt et al., 2003), which introduces additional challenges into toxicology research. Inbred mouse strains are reasonably well suited for identifying whole-genome response signatures indicative of chemical exposure because much is known regarding genetic lineage and derivation for hundreds of strains, and the number and distribution of genetic polymorphisms among mouse strains is equal to or exceeds that in the human population (Roberts et al., 2007). This approach has the added advantage of “repeat testing” in genetically identical individuals within a given strain, yielding important information regarding reproducibility of the response.
The largest group of genes identified in this study as significantly different between individuals comprised transcripts that differ in basal levels between inbred mouse strains, despite the fact that over two thirds of all strains exhibited variable degrees of liver damage. We also observed a high degree of intrastrain variability in toxicity that has been reported by other investigators, particularly for the C57BL/6J mouse (Beyer et al., 2007) and the Sprague-Dawley rat (Clayton et al., 2006). This variability can be due to a variety of factors that often cannot be controlled by the experimenter in standard toxicity studies, including epigenetic effects and differential contributions of intestinal microflora. In summary, our data underscore the value of multistrain experiments that can avert the risk of large genotype effects in a particular strain of animals used for toxicity risk assessment to determine biomarkers of response.
National Institutes of Health (U19-ES11391, P30-ES10126, and T32-ES07126); United States Environmental Protection Agency (RD832720); the U.S. Environmental Protection Agency Science to Achieve Results Graduate Fellowship (FP91690901 to A.H.).
The authors wish to acknowledge Mrs Oksana Kosyk, Mrs Blair Bradford, Dr Pierre Bushel, Dr Richard Paules, and Dr Gary A. Boorman for assistance with these studies. We thank Icoria, Inc., for the generous gift of mouse reference RNA. The research described in this article has not been subjected to each funding agency's peer review and policy review and, therefore, does not necessarily reflect their views and no official endorsement should be inferred.