Gene Expression Analysis Reveals Chemical-Specific Profiles

Hisham K. Hamadeh*, Pierre R. Bushel*, Supriya Jayadev{dagger}, Karla Martin*, Olimpia DiSorbo{dagger}, Stella Sieber*, Lee Bennett*, Raymond Tennant*,,1, Raymond Stoll{dagger}, J. Carl Barrett*, Kerry Blanchard{dagger}, Richard S. Paules* and Cynthia A. Afshari*,,1

* National Institute of Environmental Health Sciences, P.O. Box 12233, MD2-04, Research Triangle Park, North Carolina 27709; and {dagger} Boehringer-Ingelheim Pharmaceuticals, Inc., Ridgefield, Connecticut 06877

Received November 26, 2001; accepted January 8, 2002


    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
The application of gene expression profiling technology to examine multiple genes and signaling pathways simultaneously promises a significant advance in understanding toxic mechanisms to ultimately aid in protection of public health. Public and private efforts in the new field of toxicogenomics are focused on populating databases with gene expression profiles of compounds where toxicological and pathological endpoints are well characterized. The validity and utility of a toxicogenomics is dependent on whether gene expression profiles that correspond to different chemicals can be distinguished. The principal hypothesis underlying a toxicogenomic or pharmacogenomic strategy is that chemical-specific patterns of altered gene expression will be revealed using high-density microarray analysis of tissues from exposed organisms. Analyses of these patterns should allow classification of toxicants and provide important mechanistic insights. This report provides a verification of this hypothesis. Patterns of gene expression corresponding to liver tissue derived from chemically exposed rats revealed similarity in gene expression profiles between animals treated with different agents from a common class of compounds, peroxisome proliferators [clofibrate (ethyl-p-chlorophenoxyisobutyrate), Wyeth 14,643 ([4-chloro-6(2,3-xylidino)-2-pyrimidinylthio]acetic acid), and gemfibrozil (5-2[2,5-dimethylphenoxy]2-2-dimethylpentanoic acid)], but a very distinct gene expression profile was produced using a compound from another class, enzyme inducers (phenobarbital).

Key Words: gene expression; toxicogenomics; DNA arrays; classification; rat liver; pattern recognition.


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
Numerous approaches are used to investigate the relationship between chemical exposure, toxicity, and disease states. One approach is to study the modulation of gene expression in a biological model in response to chemical exposure. This modulation represents a signature of the cellular response to the effect of the studied compound. These possible signatures cannot be defined using classical methods where genes are investigated individually for potential association to chemical exposure. This is because the most highly characterized chemical-responsive genes, such as genes encoding proteins or enzymes that regulate metabolism, tend to be frequently modulated by many compounds, and therefore do not provide a solid footing for providing specificity for distinguishing multiple classes. Given the universe of compounds available, signatures may only be attained using a higher number of variables (i.e., number of genes), where the collective state (expression) of these genes would define the profile associated with exposure to that compound. The field of toxicogenomics, through the use of DNA microarrays, has the potential to advance our understanding of how multiple genes are involved in responses of biological models to chemical exposure (Burchiel 2001; Fielden and Zacharewski 2001Go; Hamadeh et al., 2001Go; Nuwaysir et al., 1999Go; Thomas et al., 2001Go; Waring et al., 2001aGo,bGo). Instead of monitoring the expression of a few genes in response to chemical exposure, DNA microarrays enable the study of levels of expression of thousands of genes at the mRNA level. The concerted expression pattern across those genes constitutes the expression profile of a compound at a certain dose and time.

Structurally unrelated compounds may belong to the same class of chemicals because of similarity in the pharmacological or toxicological endpoints they elicit. For example, at doses of diethylhexylphthalate (DEHP) and Wyeth 14,643 that produce similar levels of peroxisome proliferation in rat liver, Wyeth 14,643 produces an earlier and much greater liver tumor response than does DEHP (Biegel et al., 1992Go; Melnick et al., 1987Go).

In this study, we tested the hypothesis that structurally unrelated compounds from the same chemical class produce similar, yet distinguishable, gene expression profiles. We also hypothesized that intraclass profiles are more similar to each other than to profiles corresponding to agents from different chemical classes. In order to test whether specific patterns of gene expression can be defined for a class of compounds and whether distinguishable patterns can be discerned within that class, we studied the expression profiles of 3 well-studied agents belonging to the peroxisome proliferator class of compounds [clofibrate (ethyl-p-chlorophenoxyisobutyrate), Wyeth 14,643 (4-chloro-6[2,3-xylidino)-2-pyrimidinylthio]acetic acid) and gemfibrozil (5–2[2,5-dimethylphenoxy]2-2-dimethylpentanoic acid)]. We also studied the expression profile of a well-studied enzyme inducer, phenobarbital, in order to determine whether a distinction could be made between it and the peroxisome proliferators. Microarray analyses were performed using liver RNA derived from chemically exposed Sprague-Dawley rats at multiple time points of exposure.

This work highlights several important points for the utility of toxicogenomics studies. First, the data confirm that compound classification based on gene expression profiles is feasible. In addition, the data illustrate the differences in the gene expression elicited by chemicals that may be related in many aspects but differ with respect to toxicological effects. Further investigation of these differences might offer an explanation of the dissimilarity in adverse effects associated with various peroxisome proliferators. In addition, our report also addresses the influence of the time of exposure on gene expression by highlighting transient and delayed gene expression events in response to the 2 classes of compounds studied. As toxicogenomics databases evolve, these distinctions will be important for understanding mechanisms and developing signatures of toxicity or adaptation. Finally, the data in this paper provides important information on gene expression changes, including time-independent changes that may be used to develop signatures of the compound classes of peroxisome proliferators. Underlying the analyses of these signature genes is the potential to develop hypotheses about the potential mechanism of action of these agents.


    MATERIALS AND METHODS
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
Animal treatment and sample collection.
Male Sprague-Dawley VAF+ albino rats (CRL:CD(SD) BR; Charles River, Kingston, NY) approximately 5–7 weeks old were maintained on certified rodent chow (PMI Feeds, Inc., Brentwood, MO) ad libitum in individual stainless steel wire bottom cages suspended on racks. The animals were kept under controlled lighting (12-h light-dark cycle), temperature (72° ± 5°F), and humidity (50 ± 20%) and were acclimated to this environment for 4–7 days prior to the start of the study. Healthy rats were randomly assigned to dose groups (3 rats/group) by a computerized method. For the 2-week studies, a 1-week pretest involved dosing of all animals via oral gavage with 10 ml/kg vehicle. The 24-h studies did not have any pretest period. Clofibrate (CAS # 637–07–0), gemfibrozil (CAS # 25812–30–0), and phenobarbital (CAS # 57–30–7) were obtained from Sigma (St. Louis, MO); Wyeth 14,643 (CAS # 50892–23–4) was obtained from ChemSyn Laboratories (Lenexa, KS). Dosing suspensions of all compounds were prepared using a high speed homogenizer, and all dose suspensions were continuously stirred until completion of dosing. The peroxisome proliferators (clofibrate, gemfibrozil, and Wyeth 14,643) were prepared using 1% carboxymethylcellulose/0.2% Tween 80 as vehicle and phenobarbital was prepared using water as vehicle. Drug concentration and identity were verified via HPLC, as per United States Pharmacopeia (USP) methods for clofibrate, gemfibrozil, and phenobarbital. For Wyeth 14,643, a Waters 600E HPLC was equipped with a variable wavelength detector set at 235 nm, and a Symmetry Cs 4.6 x 250 mm column. Equal portions (approximately 10 ml) of standard and test solutions were injected separately into the column. Parameters of the run were as follows: mobile phase, 0.0738 M sodium acetate:acetonitrile (55:45 v/v); flow rate, 1 ml/min; column temperature, ambient; and run time, 6 min. The HPLC retention time of Wyeth 14,643 under these parameters was 3.4 ± 0.1 min. Body weights and food consumption were measured weekly. Based on the most recently recorded body weights, the volume of drug administered was adjusted. Animals were observed 2–3 times daily for signs of overt toxicity. Experiments were performed according to the guidelines established in the NIH Guide for the Care and Use of Laboratory Animals. At the end of the drug phase of the study, each animal was fasted overnight before necropsy. Animals were taken to a deep plain of anesthesia with CO2, sacrificed by axillary vessel incision; exsanguination and necropsy immediately followed. A cross section of the left lateral lobe of the liver was collected in 10% neutral buffered formalin for histopathology. The remaining portions of liver were collected in RNase-free tubes and snap frozen in liquid nitrogen. Frozen tissues were stored at –70°C until processed for RNA extraction. A control sample was generated by pooling livers of 9 vehicle-treated rats.

Histopathological analysis.
The liver tissues collected in formalin at necropsy were processed, embedded in paraffin, sectioned at 5 microns, and stained with hematoxylin and eosin (H&E). Histopathologic examinations of the liver sections were conducted by a pathologist and peer-reviewed.

RNA isolation.
Total RNA was isolated using QIAGEN (Qiagen, Valencia, CA) RNeasy kits. Liver sections of 130–250 mg were used for midipreps and liver sections of approximately 800 mg were used for maxipreps. Homogenization buffer was added to frozen liver sections, and the tissue was immediately homogenized on ice (tissue did not thaw prior to homogenization) using a Cyclone homogenizer equipped with a rotor/stator shaft (VirTis Company, Gardiner, NY). Homogenates were processed as per the standard QIAGEN 3/99 protocol. Final product yielded 260 nm/280 nm ratios of 1.6–2.0, purity was confirmed via gels, and concentration was determined based on 260 nm absorbances.

cDNA microarray hybridization and analysis.
A cDNA NIEHS Rat Chip, v1.0, developed in-house at NIEHS, was used for gene expression profiling experiments. A complete listing of the genes on this chip is available at the following Web site: http://dir.niehs.nih.gov/microarray/chips.htm. cDNA microarray chips were prepared according to DeRisi et al., 1996. The spotted cDNAs were derived from a collection of sequence-verified clones that covered the 3` end of the gene and ranged in size from 500 to 2000 base pairs (Research Genetics). M13 primers were used to amplify insert cDNAs from purified plasmid DNA in a 100 µl polymerase chain reaction (PCR) mixture. A sample of the PCR products (10 µl) was separated on 2% agarose gels to ensure quality of the amplifications. The remaining PCR products were purified by ethanol precipitation, resuspended in ArrayIt buffer (Telechem, San Jose CA) and spotted onto poly-L-lysine-coated glass slides using a modified, robotic DNA arrayer (Beecher Instruments, Bethesda MD).

For microarray hybridizations, each total RNA sample (35–75 µg) was labeled with Cyanine 3 (Cy3)-or Cyanine 5 (Cy5)-conjugated dUTP (Amersham, Piscataway, NJ) by a reverse transcription reaction using reverse transcriptase, SuperScript (Invitrogen, Carlsbad, California), and the primer, Oligo dT (Amersham, Piscataway, NJ). Control samples were labeled with Cy3 while samples derived from chemically exposed animals were labeled with Cy5. The fluorescently labeled cDNAs were mixed and hybridized simultaneously to the cDNA microarray chip. Each RNA pair was labeled and hybridized independently in triplicate to a total of 3 arrays. The cDNA chips were scanned with an Axon Scanner (Axon Instruments, Foster City CA) using independent laser excitation of the 2 fluors at 532 and 635 nm wavelengths for the Cy3 and Cy5 labels, respectively.

The raw pixel intensity images were analyzed using the ArraySuite v1.3 extensions of the IPLab image processing software package (Scanalytics, Fairfax, VA). This program uses methods that were developed and previously described by Chen et al. (1997) to locate targets on the array, measure local background for each target, and subtract it from the target intensity value, and to identify differentially expressed genes using a probability-based method. After pixel intensity determination and background subtraction, the ratio of the intensity of the treated cells to the intensity of the control was calculated. We have previously determined that significant autofluorescence of the gene features on the array, attributed to spotting solution, occurs at high scanning power (Tucker et al., unpublished). We measured the pixel intensity level of "blank" spots comprised of spotting solution. The data was then filtered to provide a cut off at the intensity level just above the blank measurement values in order to remove from further analyses those genes having one or more intensity values in the background range. The ratio intensity data from all of the 1700 spots printed on the NIEHS Rat Chip v1.0 was used to fit a probability distribution to the ratio intensity values and estimate the normalization constants (m and c) that this distribution provides. The constant m, which provides a measure of the intensity gain between the two channels, was approximately equal to 1 for all arrays, indicating that the channels were approximately balanced. For each array, the ratio intensity values were normalized to account for the imbalance between the 2 fluorescent dyes by multiplying the ratio intensity value by m. A probability distribution was fit to the data and used to calculate a 95% confidence interval for the ratio intensity values. Genes having normalized ratio intensity values outside of this interval were considered significantly differentially expressed.

For each of the 3 replicate arrays for each sample, lists of differentially expressed genes at the 95% confidence level were created and deposited into the NIEHS MAPS database (Bushel et al., 2001Go). For each time point and each animal, a query of the database yielded a list of genes that were differentially expressed in at least 2 of the 3 replicate hybridizations. A calculation using the binomial probability distribution indicated that the probability of a single gene appearing on this list when there was no real differential expression is approximately 0.0006. Hierarchical cluster analysis was carried out with the Cluster/TreeView package (Eisen et al., 1998Go). The entire data set is available at http://dir.niehs.nih.gov/microarray/datasets.

Real-time quantitative PCR.
RNA samples representing single animals treated with a peroxisome proliferator or phenobarbital for 24 h or 2 weeks (1852 [clofibrate, 24 hr], 1868 [Wyeth 14,643, 24 hr], 1878 [gemfibrozil, 24 hr], 1890 [phenobarbital, 24 h], 888 [clofibrate, 2 weeks], 898 [Wyeth 14,643, 2 weeks], 912 [gemfibrozil, 2 weeks], and 926 [phenobarbital, 2 weeks]) were used to validate the expression profile of 10 genes obtained using cDNA microarray data [AA818412 p450 2B2; AA996791 carnitine palmitoyl transferase 1; AI111901 tripeptidylpeptidase II; AA923966 Aflatoxin aldehyde reductase; AA957359 p55cdc; AA957519 stathmin cytosolic phosphoprotein p19; AA965078 enoyl CoA isomerase; AA818188 ketoacyl thiolase; AA963928 Ah receptor; AI070587 carboxylesterase precursor].

The primers for the aforementioned genes were designed using Primer Express® software (Applied Biosystems, Foster City, CA) and custom made (Research Genetics, Huntsville, AL). Primers that resulted in a single product which could be visualized on a 2% agarose gel were as follows: p450 2B2 [forward primer AGTGCATCACAGCCAACATCA, reverse primer GAGGGAAAAGGTCCGGTAGAA]; carboxylesterase precursor [forward primer AGTACTGGGCCAATTTTGCAA, reverse primer TGGGTGTCCAACTGCAGGTA]; Ah receptor [forward primer CATCCTGGAAATTCGAACCAA, reverse primer TGCAAGAAGCCGGAAAACT]; carnitine palmitoyl transferase 1 [forward primer CGGTTCAAGAATGGCATCATC, reverse primer ATCACACCCACCACCACGATA]; ketoacyl thiolase [forward primer ACGTGAGTGGAGGTGCCATAG, reverse primer CTCGACGCCTTAACTCGTGAAC]; stathmin p19 [forward primer CACAATCCACTGGCAAGGAA, reverse primer TGCCATGTTGGACAGAAGACA].

Real-time PCR targeting the message corresponding to these 10 genes was performed using the ABI prism 7700 Sequence Detection System (Applied Biosystems, Foster City, CA) according to the manufacturer's instructions. The SYBR® Green I labeling kit (Applied Biosystems, Foster City, CA), was used to detect double-stranded DNA generated during PCR amplification, used according to the manufacturer's instructions. Reverse transcription and PCR reactions were performed at the same time in a 50 µl reaction containing 4 mM MgCl2, 0.8 mM of each dNTP, 100 ng total RNA, 0.4 µM reverse primer and 0.4 µM forward primer, 0.4 units/µl RNasin, 0.025 units/µl AmpliTaq Gold DNA polymerase (Roche, Basel, Switzerland) and 0.25 units/µl MulV Reverse Transcriptase (Roche, Basel, Switzerland). Amplification reactions were carried out using the following temperature profile: 48°C, 30 min; 95°C, 10 min; 95°C, 15 s; 60°C, 1 min) for 40 cycles. Fluorescence emission was detected for each PCR cycle and the threshold cycle (CT) values were determined. The CT value was defined as the actual PCR cycle when the fluorescence signal increased above the background threshold. Induction or repression of a gene in a treated sample relative to control was calculated as follows: Fold increase/decrease = e – (CT(exposed) – CT(control)). Values were reported as an average of triplicate analyses.


    RESULTS
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
Histopathology.
To investigate whether different chemical agents produce distinguishable gene expression profiles in biological systems, compounds belonging to 2 classes of rodent liver toxicants were chosen for study (IARC 1977Go,1987Go). Sprague-Dawley rats were exposed to 1 of 3 peroxisome proliferators (clofibrate, Wyeth 14,643, or gemfibrozil), or to an enzyme inducer (phenobarbital), as described in Materials and Methods. Animals were dosed orally via gavage with 250 mg/kg/day of clofibrate, 250 mg/kg/day of Wyeth 14,643, 100 mg/kg/day of gemfibrozil, or 120 mg/kg/day of phenobarbital in a volume of 10 ml/kg. The administered doses targeted the maximum tolerated dose (MTD) for each of the chemicals (Butterworth et al., 1995Go). Histopathological examinations of liver sections were conducted as described in the methods section to score the gross tissue and organ effect of the administered doses at each time. No drug-related microscopic observations were apparent in animals sacrificed 24 h after a single treatment with any of the compounds, whereas drug-related microscopic hepatocellular hypertrophy was observed in livers of all the 2-week treated animals. For the 3 peroxisome proliferators, hypertrophic hepatocytes were characterized by large cells with abundant microvesiculated eosinophilic cytoplasm. Large cells with abundant pale-staining eosinophilic cytoplasm with basophilic stippling characterized hypertrophic hepatocytes, in animals treated with phenobarbital.

Gene expression.
In order to determine gene expression changes associated with chemical exposure, liver mRNA was collected 24 h following a single exposure or after 2 weeks of daily exposures to the compound. Competitive hybridizations of fluorescently labeled cDNA (Duggan et al., 1999Go) derived from control vs. treated livers were used to measure relative abundance of each mRNA on the custom NIEHS cDNA microarray Rat Chip v1.0, which contained ~1700 sequence-verified rat genes. We conducted statistical analyses of the microarray data and determined significantly changed genes using a ratio distribution model at the 95% confidence level, as mentioned in the Materials and Methods section. We were able to reduce the probability of false positives in our data set by performing triplicate hybridizations on each of at least 3 independent biological samples and by including only genes exhibiting a binomial distribution probability <= 0.0006 (Bushel et al., 2001Go). These genes were utilized for further higher-level comparative analyses (e.g., clustering).

Exploratory interpretation.
The results obtained from the collective microarray analysis of the peroxisome proliferator-treated rat livers revealed significant gene expression changes in approximately 25% of the genes on the rat chip and elucidated interesting molecular changes and pathway relationships associated with peroxisome proliferator exposure (Table 1Go). These pathways include stimulation of triglyceride hydrolysis, fatty acid uptake and conversion to acyl CoA derivatives, and stimulation of the ß-oxidation pathway. Observation of alteration of these pathways corroborates past data (Amacher et al., 1997Go; Schoonjans et al., 1996Go) and serves as a validation of the use of microarrays to rapidly interrogate effector pathways of toxicants.


View this table:
[in this window]
[in a new window]
 
TABLE 1 Peroxisome Proliferator Effects on a Sampling of Genes
 
The mechanism of action of phenobarbital, a compound that has been studied for over 40 years, is only partially understood. Similar to the peroxisome proliferator data, our phenobarbital microarray data corroborates previously described metabolic, pharmacologic, and toxicologic effects of phenobarbital and offers new clues that might be further investigated to better define its mechanism of action (Table 2Go). For example, we observed the upregulation of previously reported cytochrome P450 genes, such as CYP 2B2, 2C6, 3A9, and other genes such as epoxide hydrolase, diaphorase, and several GSTs (Furukawa et al., 1985Go; Griffin et al., 1984Go; Tavoloni et al., 1983Go; Whysner et al., 1996Go) as well as the induction of several novel genes, such as carboxylesterase precursor. In addition, our novel observation of the downregulation of carnitine palmitoyl transferase 1 (CPT 1) in phenobarbital-treated rats may explain the significant reduction (30–60%) of 4 carnitine constituents (total and free carnitine and short-and long-chain fatty acid carnitine esters) observed in serum from 471 patients treated for convulsions with phenobarbital (Hug et al., 1991Go). The combination of decreased CPT 1 levels and downregulation of Acyl CoA synthetase, a gene involved in the catalysis of fatty acid esterification, suggests an inhibition of fatty acid peroxidation. Major metabolic pathways such as gluconeogenesis, glycolysis, ß-oxidation, and fatty acid peroxidation were decreased by phenobarbital treatment, while fatty acid synthesis was stimulated (Argaud et al., 1991Go; Thurman and Marazzo, 1975Go). Although some of these gene inductions were previously reported, there has been little effort to integrate the observations in the context of phenobarbital's mechanism of action. Our results offer a more comprehensive overview of molecular responses to toxicant exposure by revealing the coordinate expression of multiple genes in homeostasis and metabolic pathways. The agreement of the expression profiles for phenobarbital and peroxisome proliferators with past, traditional studies lends confidence in the use of these gene expression profiles in further pattern recognition applications.


View this table:
[in this window]
[in a new window]
 
TABLE 2 Phenobarbital Effects on a Sampling of Genes
 
Gene expression validation.
Our gene expression data is validated in 3 ways. The first is by replicate analysis. We used 3 animals for each compound and measured the gene expression for each animal on 3 chips. This approach generated 9 measurements for each gene. The combined use of a confidence interval and the binomial probability aided in eliminating genes with high biological or technical variability from further analyses (i.e., clustering). Secondly, we routinely conduct resequencing of the clones we find significantly changed. Currently, our clone set shows an accuracy of approximately 90%. Identification of clones is updated on our Web site, http://dir.niehs.nih.gov/microarray/chips.htm. Finally, we validated the expression profile of 10 genes across samples derived from individual animals exposed to peroxisome proliferators or phenobarbital at the 24-h and 2-week time points. Sample PCR products were run on a 2% agarose gel and visualized by ethidium bromide staining. A representative gel indicating the quality of the reactions is shown in Figure 1Go. Comparison of data from cDNA microarray and real-time polymerase chain reaction (RT-PCR) evaluations demonstrated a high level of correlation between the 2 approaches (Table 3Go), where the induction or repression of each gene was confirmed across multiple samples. Microarray measurements are typically only semiquantitative, with compression of values occurring at high-fold changes. The RT-PCR measurements are likely to provide better quantitation for genes such as p4502B2 (phenobarbital 12.42 fold on microarray, 32-fold by RT-PCR), but generally the quantitative measurements with the 2 approaches are well correlated (Table 3Go).



View larger version (63K):
[in this window]
[in a new window]
 
FIG. 1. Assessment of amplification of the RT-PCR specificity validation experiment using agarose gel electrophoresis: 60 µl of PCR reaction product corresponding to each gene were concentrated, then separated on a 2% agarose gel. Marker ladders (1 kb and 100 bp; Invitrogen, Carlsbad, CA) were run on either side of the samples. Each lane shows a representative product for the gene indicated.

 

View this table:
[in this window]
[in a new window]
 
TABLE 3 Validation of cDNA Microarray Data by Real-Time PCR
 
Pattern recognition.
A critical question in toxicogenomics is whether gene expression information may be used to reveal chemical-specific signature patterns. We used several computational analyses to determine whether different toxicants result in distinguishable gene expression patterns. Application of hierarchical cluster analysis (Eisen et al., 1998Go) confirmed that individual animals could be distinguished by the class of toxicants to which they were exposed (Fig. 2Go) and revealed 2 distinct nodes containing animals treated with either of the 2 classes of chemicals. Across experiments, one node represents the 3 animals that were exposed to phenobarbital, while the other node includes animals treated with Wyeth 14,643, gemfibrozil, or clofibrate. Furthermore, in the latter node, individual animals were clustered in separate subnodes denoting the specific peroxisome proliferator to which they were exposed.



View larger version (49K):
[in this window]
[in a new window]
 
FIG. 2. Different class compound toxicants generate discrete gene expression profiles. Genes whose expression profile was significantly altered in any toxicant-exposed animal at the 24-h time point were clustered according to their expression levels across animals. Samples were also clustered based on the pattern of expression of the studied genes. A node is highlighted showing the expression pattern of a subset of genes representative of the overall grouping of samples derived from the livers of chemically exposed rats.

 
We used principal components analysis (PCA), which is a pattern recognition technique that represents a multivariate statistical method that is useful for reducing multidimensional data (such as high-density gene expression data) down to 2 or 3 dimensions that can be readily comprehended. The principal components were new variables created from linear combinations of the starting variables (genes), where each principal component is orthogonal or not correlated with all others. The first principal component contained the largest part of the variance of the data set (37.5%) with the subsequent principal components containing correspondingly smaller amounts of variance (18.7 and 12.1% for 2nd and 3rd, respectively). This analysis allowed the visualization, in 3-dimensional space for simplicity, of the discrimination between the gene expression responses elicited by these 2 classes of compounds (Fig. 3Go). PCA demonstrated close proximity in the gene expression pattern between clofibrate-, Wyeth 14,643-, and gemfibrozil-exposed animals, but indicated a distinct partition (separation) between these compounds and phenobarbital-exposed animals.



View larger version (12K):
[in this window]
[in a new window]
 
FIG. 3. The Partek Pro 2000 software package was used for visual PCA of the data for genes that were altered in a statistically significant manner with any of the treatments used. The first principal component for this data set accounted for 37.4% of the variation present, the second component for an additional 18.7%, and the third for 12.1%. Each colored point represents data from an individual animal treated with the respective agent for 24 h.

 
Finally, pairwise comparisons (Fig. 4Go) of gene expression profiles of individual animals exposed to chemicals showed confirmation of the discerning potential of microarray data. Comparison of gene expression profiles of two different animals exposed to the same compound resulted in a relatively high correlation (e.g., r > +0.8 for Wyeth 14,643 vs. Wyeth 14,643) as compared to animals exposed to different compounds belonging to the same class (e.g., r ~ +0.6 for Wyeth 14,643 vs. gemfibrozil or clofibrate). However, comparisons of animals treated with toxicants belonging to different classes of compounds resulted in a relatively low correlation (e.g., r < + 0.4 for Wyeth 14,643 vs. phenobarbital). In summary, the data in Figures 2, 3, and 4GoGoGo indicate that through multiple approaches and bioinformatics tools, it is possible to discriminate gene expression profiles generated as a response to distinct liver toxicants. Proof of this concept aids in validation of the future potential of a predictive toxicogenomic strategy.



View larger version (34K):
[in this window]
[in a new window]
 
FIG. 4. The ratios of exposed to control animals corresponding to transcript levels of genes whose expression was significantly altered in any toxicant-exposed animal at the 24-h time point, were compared via set pairwise correlation analysis using Spotfire software (Spotfire, Inc., Cambridge, MA). Comparison of rats treated with the sample compound showed the highest correlation (Wyeth 14,643-treated animals). Comparison of Wyeth 14,643 with either clofibrate-or gemfibrozil-treated animals yielded appreciable correlation, which was in agreement with the fact that both of those compounds belonged to the same peroxisome proliferator class of toxicants. The correlation dropped sharply when Wyeth 14,643-and phenobarbital-treated animals were compared, and indicated poor correlation.

 
Time dependency of gene expression.
One potential application of microarrays in toxicology is their use in predicting toxicity of undefined chemicals prior to the appearance of pathological or disease outcomes. Gene expression profiles from animals treated with a compound of unknown toxicity may be compared with a database of DNA microarray-generated gene expression data corresponding to known compounds. However, toxicological effects are confounded by time and dose dependency of lesions that lead to differences in gene expression signatures. Therefore, we evaluated whether time-independent gene expression responses that represent signatures of chemical-specific alterations occur. Application of clustering methods to the data from phenobarbital and the peroxisome proliferators at 2 time points allowed the identification of time-independent regulation (up or down) of genes in response to compound treatment (Fig. 5Go). Those genes that are regulated in the same fashion upon compound exposure at multiple time points may potentially serve as reliable biomarkers of effect when establishing time-independent gene expression profiles.



View larger version (49K):
[in this window]
[in a new window]
 
FIG. 5. Two-dimensional, hierarchical clustering analysis of genes that were altered in a statistically significant manner at 95% confidence level in at least 2 of the replicate hybridizations performed on each sample. Data from 9 hybridizations, representing 3 replicates of 3 independent biological samples derived from rats treated with the same compound at each time point, were averaged, and those values were used for clustering analysis. Clustering analysis was performed across genes but not animals. Genes are represented on the vertical axis while animals are on the horizontal axis. (A) Highlighted nodes from a hierarchical clustering tree, showing a suite of genes whose transcripts were increased in phenobarbital-exposed animals over controls in a time-independent fashion. (B) Genes modulated by gemfibrozil treatment in a similar fashion at 24 h and 2 weeks of exposure. (C) Nodes showing different genes that were modulated by peroxisome proliferators collectively in a time-independent fashion. Red, gene induction; green, gene repression in the treated samples, relative to control.

 
RNA based gene expression analysis may be regarded as a snapshot of molecular occurrences in time/space. Analysis of expression at one time point can be misleading due to transient and delayed gene expression events. We studied transient and delayed alterations in gene expression in response to the daily exposure to the chemicals used in this study. Figure 6Go illustrates the value of studying multiple time points, and the cluster indicates genes that are altered at the delayed time point. These genes may be further investigated or corroborated in future studies for their association with chronic toxicity or adaptation (Fig. 6Go).



View larger version (48K):
[in this window]
[in a new window]
 
FIG. 6. Illustration of transient versus delayed responses in gene expression. Nodes from 4 hierarchical clustering trees corresponding to the chemicals used in this study were highlighted to show examples of (A) transiently altered transcripts or (B) genes that required a delayed period of time for up-or downregulation, implicating those genes as possible biomarkers of toxicant-associated delayed toxicity or adaptation to exposure. Red, gene induction; green, gene repression in the treated samples, relative to control.

 

    DISCUSSION
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
The goal of this study was to determine whether generation of chemically associated gene expression profiles, using microarray technology, would permit classification of compound associated signatures. We generated gene expression profiles corresponding to 3 peroxisome proliferators, clofibrate, Wyeth 14,643, and gemfibrozil, and an enzyme inducer, phenobarbital, and demonstrated that gene expression profiling is indeed a powerful tool for distinguishing gene expression generated by structurally unrelated toxicants in an in vivo model. Those distinctions were made even when compounds shared some endpoints such as peroxisome proliferation. A greater similarity was found among peroxisome proliferators than among the peroxisome proliferators and phenobarbital. Clustering of genes that were significantly affected by the 24-h exposures, where no histopathological abnormalities were detected, demonstrated that gene expression profiles might be successfully used for compound classification.

Whether animals should be grouped together as a pool or examined individually represents one issue in the design of toxicogenomics studies. Some investigators advocate pooling for the costly microarray analysis and using individual animals for the follow-up verification steps. However, pooling may cause misinterpretation of data if one animal shows a remarkably distinct response, or lack of response. In this study we analyzed individual chemically exposed animals against a pool of control animals. The generation of gene expression profiles corresponding to 3 animals exposed to the same compound, as opposed to pooling, allowed for the detection of interanimal differences. This facilitated the testing of the robustness of DNA microarray technology and pattern recognition algorithms to generate distinguishable gene expression profiles, despite the existence of differences in response across similarly treated animals. Although we have observed similarities as high as about 90% between animals exposed to the same agents, this similarity showed chemical dependence, reaching about 80% with some of the compounds. The best interanimal correlation was observed among Wyeth 14,643-treated animals, where we saw the greatest number of statistically significant gene expression modulations, suggesting that compound potency may be positively related to decreased variation in gene expression response.

There is high concordance of the expression changes found in the microarray analyses in phenobarbital-exposed rats with the results obtained by scientists using other methodologies over many years of study. This concordance is illustrated in a display for gene expression changes in the form of an "effector pathway map" (EPM) for chemical action. The phenobarbital gene expression profile data set was mapped into previously defined cellular response pathways to demonstrate the informational power of this type of data presentation (Fig. 7Go). Organization of data in this integrated diagram facilitates clear visualization of regulatory and molecular events occurring as a response to chemical exposure and visualization of pathways that were affected by phenobarbital. At a glance, this diagram illustrates how the gene expression profile data has corroborated past findings (Fig. 7Go, yellow targets) and may contribute new mechanistic insights (Fig. 7Go, blue targets).



View larger version (44K):
[in this window]
[in a new window]
 
FIG. 7. Genes modulated with phenobarbital exposure, displayed as a map of effector pathways of toxicant (MEPT), a schematic diagram of phenobarbital-modulated gene expression events. Enzymes enclosed in yellow boxes were altered in this study and were previously reported in the literature to be modulated by phenobarbital treatment, while those in blue were not previously associated with phenobarbital exposure and are novel observations. The colored circles indicate the type of modulation. The circle before the slash in each grouping denotes the 24-h time point; the circle after the slash, the 2-week time point. Red, upregulation; green, downregulation; black, no change.

 
Chronic cell proliferation is proposed as a major mechanism for tumor promotion by phenobarbital (Barbason et al., 1983Go; Feldman et al., 1981Go; Whysner et al., 1996Go). Histopathological analysis indicated liver enlargement at 2 weeks in phenobarbital-treated animals. Several studies have demonstrated that phenobarbital induced increases in DNA synthesis in rats and mice of various strains (Busser and Lutz, 1987Go). Our data corroborates phenobarbital-induced increase in cell proliferation. Transcript levels for DNA polymerase b, AIRC-SAICAR synthase, cyclin B1, and ARF-like factor (ARL5) were increased at the 24-h and 2-week time points, relative to controls, by phenobarbital exposure, indicating increased DNA synthesis and liver proliferation. This observation is further supported by the modulation in cytoskeleton rearrangement-related genes, such as rho-associated protein kinase, which is essential for the rearrangement of actin cytoskeleton (Ohashi et al., 2000Go; Watanabe et al., 1999Go), suggesting a role for these alterations in the hypertrophy observed in the livers of the exposed animals. We also observed the upregulation of phosphoprotein stathmin (p19) at the 2 time points in response to phenobarbital exposure. Phosphoprotein stathmin (p19), which is abundant in neuroendocrine tumor cells, showed a 15-fold greater abundance in newborn than in adult brain and its levels increased after two-thirds partial rat hepatectomy (Koppel et al., 1993Go), suggesting a role in regeneration and growth of various tissues. These data support the observation that phenobarbital produces liver enlargement due to proliferation of liver cells and offers new insights on its molecular basis.

Analysis of transient, delayed, and time-independent alterations suggested a relationship between the pattern of expression and gene function. The majority of genes that were expressed in a time-independent fashion in response to phenobarbital or peroxisome proliferator treatment corresponded to enzymes that function in compound metabolism (Fig. 5AGo; p450 2B2, epoxide hydrolase 1, GST, aflatoxin B1 aldehyde reductase) or cell biochemical processes (Figs. 5A–5C Go;Acyl CoA oxidase, Carnitine palmitoyltransferase 1, histidine decarboxylase, stearyl-CoA desaturase). This makes sense when one considers that animals were treated with the studied chemicals on a daily basis, furnishing recurring surges in blood levels of the compounds in the exposed animals, and thus affecting compound metabolism genes or downstream effects.

Metabolism-related genes were notably absent from the transiently altered transcripts, the majority of which represented signaling related genes (Fig. 6AGo). This is consistent with the transient nature of the response and these genes probably constituted an initial response in the liver upon exposure to the specific toxicants for the first time. Alterations in gene expression that were present at 2 weeks of exposure but not at the 24-h time point (Fig. 6BGo) constituted delayed responses to toxicant exposure. These responses could be associated with adaptation events or with the relation to the histopathological observation of hypertrophy noted in all animals treated by chemicals for 2 weeks. These changes could also be due to overt toxicity that may be manifested in gene expression responses but not necessarily detected by histopathological examination until a much longer period of exposure.

We have successfully generated gene expression profiles for 3 peroxisome proliferators and an enzyme inducer, and we were able to show, through the application of pattern recognition algorithms and computational analyses, that these patterns were distinct. We demonstrated that chemicals from the same class of compounds give rise to discernable gene expression profiles that bear more similarity to each other than to patterns corresponding to exposure to compounds from a different class. Systematic development of an expanded database for gene expression responses to drugs and environmental pollutants may yield compound-specific signature patterns that would also provide insights into affected regulatory and proliferative and repair/adaptive pathways. We demonstrated the validity of our expression data by corroborating published reports on the chemicals that we utilized. In addition, our data revealed gene expression that has not been previously associated with the compounds we used and suggest that such results will provide valuable targets for further investigations of the mechanism of action of chemical hazards.


    ACKNOWLEDGMENTS
 
We thank Julie Foley, Jeffrey Tucker, Kelli Reynolds, Jeffrey Reece, Jennifer Collins, and Dr. Johnathan Ladapo for assistance in analyzing and compiling the data from this study. In addition, we acknowledge and thank Drs. Jeff Trent, Mike Bittner, and Pat Hurban, and Donald Cox, for advice in development of our microarray center. We also thank Nigel Walker and Chris Miller for their help and guidance in the RT-PCR analysis. Finally, we thank Drs. Thomas Kunkel, Robert Maronpot, Samuel Wilson, Paul Nettesheim, and Masahiko Negishi for helpful comments on the preparation of this manuscript.


    NOTES
 
1 To whom correspondence should be addressed. Fax: (919) 316-4535. E-mail: afshari{at}niehs.nih.gov. Back


    REFERENCES
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS
 DISCUSSION
 REFERENCES
 
Amacher, D. E., Beck, R., Schomaker, S. J., and Kenny, C. V. (1997). Hepatic microsomal enzyme induction, ß-oxidation, and cell proliferation following administration of clofibrate, gemfibrozil, or bezafibrate in the CD rat. Toxicol. Appl. Pharmacol. 142, 143–150.[ISI][Medline]

Argaud, D., Halimi, S., Catelloni, F., and Leverve, X. M. (1991). Inhibition of gluconeogenesis in isolated rat hepatocytes after chronic treatment with phenobarbital. Biochem. J. 280, 663–669.[ISI][Medline]

Barbason, H., Rassenfosse, C., and Betz, E. H. (1983). Promotion mechanism of phenobarbital and partial hepatectomy in DENA hepatocarcinogenesis cell-kinetics effect. Br. J. Cancer 47, 517–525.[ISI][Medline]

Biegel, L. B., Hurtt, M. E., Frame, S. R., Applegate, M., O`Connor, J. C., and Cook, J. C. (1992). Comparison of the effects of Wyeth 14,643 in Crl:CD BR and Fisher-344 rats. Fundam. Appl. Toxicol. 19, 590–597.[ISI][Medline]

Burchiel, S. K., Knall, C. M. Davis, J. W., III, Paules, R. S. Boggs, S. E. Afshari, C. A. (2001). Analysis of genetic and epigenetic mechanisms of toxicity: Potential roles of toxicogenomics and proteomics in toxicology. Toxicol. Sci. 59, 193–195.[Abstract/Free Full Text]

Bushel, P. R., Hamadeh, H., Bennett, L., Sieber, S., Martin, K., Nuwaysir, E. F., Johnson, K., Reynolds, K., Paules, R. S., and Afshari, C. A. (2001). MAPS: A microarray project system for gene expression experiment information and data validation. Bioinformatics 17, 564–565.[Abstract/Free Full Text]

Busser, M. T., and Lutz, W. K. (1987). Stimulation of DNA synthesis in rat and mouse liver by various tumor promoters. Carcinogenesis 8, 1433–1437.[Abstract]

Butterworth, B. E., Conolly, R. B., and Morgan, K. T. (1995). A strategy for establishing mode of action of chemical carcinogens as a guide for approaches to risk assessments. Cancer Lett. 93, 129–146.[ISI][Medline]

Chen, Y., Dougherty, E, R., and Bittner, M. L. (1997). Ratio-based decisions and the quantitative analysis of cDNA microarray images. J. Biomed. Optics 2, 364–374.

DeRisi, J., Penland, L., Brown, P. O., Bittner, M. L., Meltzer, P. S., Ray, M., Chen, Y., Su, Y. A., and Trent, J. M. (1996). Use of a cDNA microarray to analyse gene expression patterns in human cancer. Nat. Genet. 14, 457–460.[ISI][Medline]

Duggan, D. J., Bittner, M., Chen, Y., Meltzer, P., and Trent, J. M. (1999). Expression profiling using cDNA microarrays. Nat. Genet. 21(Suppl. 1), 10–14.[ISI][Medline]

Eisen, M. B., Spellman, P. T., Brown, P. O., and Botstein, D. (1998). Cluster analysis and display of genome-wide expression patterns. Proc. Natl. Acad. Sci. U.S.A 95, 14863–14868.[Abstract/Free Full Text]

Feldman, D., Swarm, R. L., and Becker, J. (1981). Ultrastructural study of rat liver and liver neoplasms after long-term treatment with phenobarbital. Cancer Res. 41, 2151–2162.[Abstract]

Fielden, M. R., and Zacharewski, T. R. (2001). Challenges and limitations of gene expression profiling in mechanistic and predictive toxicology. Toxicol. Sci. 60, 6–10.[Abstract/Free Full Text]

Furukawa, K., Numoto, S., Furuya, K., Furukawa, N. T., and Williams, G. M. (1985). Effects of the hepatocarcinogen nafenopin, a peroxisome proliferator, on the activities of rat liver glutathione-requiring enzymes and catalase in comparison to the action of phenobarbital. Cancer Res. 45, 5011–5019.[Abstract]

Griffin, M. J., Kirsten, E., Carubelli, R., Palakodety, R. B., McLick, J., and Kun, E. (1984). The in vivo effect of benzamide and phenobarbital on liver enzymes: Poly(ADP-ribose) polymerase, cytochrome P-450, styrene oxide hydrolase, cholesterol oxide hydrolase, glutathione S-transferase and UDP-glucuronyl transferase. Biochem. Biophys. Res. Commun. 122, 770–775.[ISI][Medline]

Hamadeh, H. K., Bushel, P. B., Paules, R., and Afshari, C. A. (2001). Discovery in toxicology: Mediation by gene expression array technology. J. Biochem. Mol. Toxicol. 15, 231–242.[ISI][Medline]

Hug, G., McGraw, C. A., Bates, S. R., and Landrigan, E. A. (1991). Reduction of serum carnitine concentrations during anticonvulsant therapy with phenobarbital, valproic acid, phenytoin, and carbamazepine in children. J. Pediatr. 119, 799–802.[ISI][Medline]

IARC (1977). IARC monographs on the evaluation of the carcinogenic risk of chemicals to man: Some miscellaneous pharmaceutical substances. IARC Monogr. Eval. Carcinog. Risk Chem. Man. 13, 1–255.[Medline]

IARC (1987).Overall evaluations of carcinogenicity: An updating of IARC Monographs volumes 1 to 42. IARC Monogr. Eval. Carcinog. Risks Hum. 7 (Suppl.), 1–440.

Koppel, J., Loyer, P., Maucuer, A., Rehak, P., Manceau, V., Guguen-Guillouzo, C., and Sobel, A. (1993). Induction of stathmin expression during liver regeneration. FEBS Lett. 331, 65–70.[ISI][Medline]

Melnick, R. L., Morrissey, R. E., and Tomaszewski, K. E. (1987). Studies by the National Toxicology Program on di(2-ethylhexyl)phthalate. Toxicol. Ind. Health 3, 99–118.[ISI][Medline]

Nuwaysir, E. F., Bittner, M., Trent, J., Barrett, J. C., and Afshari, C. A. (1999). Microarrays and toxicology: The advent of toxicogenomics. Mol. Carcinog. 24, 153–159.[ISI][Medline]

Ohashi, K., Nagata, K., Maekawa, M., Ishizaki, T., Narumiya, S., and Mizuno, K. (2000). Rho-associated kinase ROCK activates LIM-kinase 1 by phosphorylation at threonine 508 within the activation loop. J. Biol. Chem. 275, 3577–3582.[Abstract/Free Full Text]

Schoonjans, K., Staels, B., and Auwerx, J. (1996). Role of the peroxisome proliferator-activated receptor (PPAR) in mediating the effects of fibrates and fatty acids on gene expression. J. Lipid Res. 37, 907–925.[Abstract]

Tavoloni, N., Jones, M. J., and Berk, P. D. (1983). Dose-related effects of phenobarbital on hepatic microsomal enzymes. Proc. Soc. Exp. Biol. Med. 174, 20–27.[Abstract]

Thomas, R. S., Rank, D. R., Penn, S. G., Zastrow, G. M., Hayes, K. R., Pande, K., Glover, E., Silander, T., Craven, M. W., Reddy, J. K., Jovanovich, S. B., and Bradfield, C. A. (2001). Identification of toxicologically predictive gene sets using cDNA microarrays. Mol. Pharmacol. 60, 1189–1194.[Abstract/Free Full Text]

Thurman, R. G., and Marazzo, D. P. (1975). Mixed-function oxidation and intermediary metabolism: Metabolic interdependencies in the liver. Adv. Exp. Med. Biol. 58, 355–367.[Medline]

Waring, J. F., Ciurlionis, R., Jolly, R.A., Heindel, M., and Ulrich, R. G. (2001a). Microarray analysis of hepatotoxins in vitro reveals a correlation between gene expression profiles and mechanisms of toxicity. Toxicol. Lett. 120, 359–368.[ISI][Medline]

Waring, J. F., Jolly, R. A., Ciurlionis, R., Lum, P. Y., Praestgaard, J. T., Morfitt, D. C., Buratto, B., Roberts, C., Schadt, E., and Ulrich, R. G. (2001b). Clustering of hepatotoxins based on mechanism of toxicity using gene expression profiles. Toxicol. Appl. Pharmacol. 175, 28–42.[ISI][Medline]

Watanabe, N., Kato, T., Fujita, A., Ishizaki, T., and Narumiya, S. (1999). Cooperation between mDia1 and ROCK in Rho-induced actin reorganization. Nat. Cell Biol. 1, 136–143.[ISI][Medline]

Whysner, J., Ross, P. M., and Williams, G. M. (1996). Phenobarbital mechanistic data and risk assessment: Enzyme induction, enhanced cell proliferation, and tumor promotion. Pharmacol. Ther. 71, 153–191.[ISI][Medline]