Variation in the Hepatic Gene Expression in Individual Male Fischer Rats

Abstract

A new tool beginning to have wider application in toxicology studies is transcript profiling using microarrays. Microarrays provide an opportunity to directly compare transcript populations in the tissues of chemical-exposed and unexposed animals. While several studies have addressed variation between microarray platforms and between different laboratories, much less effort has been directed toward individual animal differences especially among control animals where RNA samples are usually pooled. Estimation of the variation in gene expression in tissues from untreated animals is essential for the recognition and interpretation of subtle changes associated with chemical exposure. In this study hepatic gene expression as well as standard toxicological parameters were evaluated in 24 rats receiving vehicle only in 2 independent experiments. Unsupervised clustering demonstrated some individual variation but supervised clustering suggested that differentially expressed genes were generally random. The level of hepatic gene expression under carefully controlled study conditions is less than 1.5-fold for most genes. The impact of individual animal variability on microarray data can be minimized through experimental design.

Introduction

High-density microarrays are powerful tools to identify gene expression changes and are beginning to impact toxicology by associating gene patterns with specific toxicants (Burczynski et al., 2000; Bulera et al., 2001; Holmes et al., 2001; Hamadeh et al., 2002; Gant et al., 2003; Waring et al., 2003). Gene expression measurements have largely focused on the liver because the liver is a major organ in toxicant metabolism (Ulrich et al., 2004) and because hepatotoxicity is the primary reason for withdrawal of new drugs from the market (Zimmerman, 1998; Gluud, 2002; Goodman, 2002; Teoh and Farrell, 2003). One of the benefits of microarrays is the ability to detect very small differences in transcript abundance, but this also means that they have exquisite sensitivity to detect experimental variation in gene expression (Causton et al., 2003).

The variation inherent in the microarray data can be roughly classified into two types: technical and biological. Technical variation is related to the collection and analysis of the microarray data and varies with different hybridization platforms (Causton et al., 2003; Goodsaid et al., 2004). There are several means for dealing with technical variability including replicate reporter elements on the array, dye reversal hybridizations (fluor flips), multiple reporters representing different regions of the same gene, and replicate arrays to create duplicate data sets (Nguyen et al., 2002; Causton et al., 2003; Frederiksen et al., 2003; Rosenzweig et al., 2004). Proper statistical methods can also substantially increase the chances that the gene expression changes detected represent true differential expression (Wolfinger et al., 2001; Black and Doerge, 2002; Cho and Lee, 2004).

A second source of variation in microarray data is biological variation (Causton et al., 2003). While several recent studies have addressed interlaboratory variation in microarray results, comparisons were between pooled samples (Baker et al., 2004; Chu et al., 2004; Ulrich et al., 2004; Waring et al., 2004). Thus, individual variation was not addressed. Characterizing individual variation is essential for recognizing changes associated with disease states (Lee et al., 2002; Whitney et al., 2003).

The ability to make direct comparisons between 2 samples on the same array for the same gene is a unique and important feature of the 2-color microarray (Churchill, 2002). The labeled samples from treated animals are compared with either a paired control (Gant et al., 2003; Smith et al., 2003) or a pooled control. (Hamadeh et al., 2002; Waring et al., 2003). Thus, both treated individuals as well as the untreated controls may contribute to the variability.

Many factors have the potential to affect gene expression in rodents including feeding, caging, handling, circadian cycle, and sequence of necropsy (Dunn and Scheuing, 1971; Pritchard et al., 2001; Causton et al., 2003). When rats are housed 3 or more per cage, the dominant rat may eat longer (Sharp et al., 2002, 2003). Gene expression in the liver is highly dependent on the time and amount of food ingested (Delzenne et al., 2001; Kato and Kimura, 2003).

A survey of the individual variation in gene expression patterns in the livers from control rats was done using oligonucleotide microarrays in 2 independent experiments. Environmental factors such as light intensity, temperature, humidity, type of feed, cages and bedding, housing, and animal handling procedures were closely controlled. Gene expression profiles were analyzed against known biological parameters such as body weight, histopathology, and clinical chemistry results. The magnitude of differential gene expression in these animals was small with few gene changes greater than 2-fold. This data provides a baseline for judging subtle effects from chemical exposure.

Materials and Methods

Animals

Male Fischer 344 rats, approximately 36 ± 3 days old, were supplied by Taconic laboratory animals (Germantown, NY) and were approximately 89 ±3 days old (experiment 1) or 91 ±3 days old (experiment 2) when receiving the control vehicle. The studies were conducted at Battelle International, Inc., Columbus, Ohio and the protocol was approved by the Battelle IACUC and followed the standards outlined in the Guide for the Care and Use of Laboratory Animals (NRC, 1996). Sentinel animals included in the study were negative for rodent pathogens (Rao et al., 1989). Rats were randomized to experiments by body weight partitioning using the PATH/TOX SYSTEM (Xybion Medical Systems Corp., Cedar Knolls, NJ) algorithm. The rats were housed 3 per cage in 22″ L × 12.5″ W × 8″ H polycarbonate cages (Lab Products, Inc., Seaford, DE) with polyester cage filters (Snow Filtration Co., Cincinnati, OH). Animal room temperature and humidity were continuously monitored and varied between 71 and 75°F and 36% to 48% relative humidity. Each experiment was divided into 2 lighting (acclimated for 2 weeks prior to study start) exposure groups. The rats in the day group had a 12-hour light period from 8 AM to 8 PM with a corresponding 12-hour dark (night) period from 8 PM to 8 AM. Rats in the light reversal group had a 12-hour light period from 8 PM to 8 AM with a corresponding 12-hour dark (night) period from 8 AM to 8 PM. The times and equivalent times for both day and light reversal groups are shown in Figure 1.

Serum melatonin concentrations were used to determine whether rats had successfully adapted to the light reversal. Sera were collected for melatonin analyses from a group of 4 rats at 4 PM and another 4 of rats at 6 PM 1 day prior to experiment 1 and repeated for experiment 2 in both the day and light reversal rats for a total of 32 rats. The mean (±SD) melatonin levels were 10 ±5 and 12 ±3 for the 2 samples taken during the day and 307 ±71 and 348 ±35 pg/ml for the 2 samples taken during the dark in the light reversal rats. The melatonin levels demonstrate that the rats had successfully adapted to the light reversal (Travlos et al., 2001).

The rats had ad libitum access to irradiated NTP-2000 wafer feed (Ziegler Brothers, Gardners, PA) during the 12-hour dark period with no food present in their cage during the 12-hour light period of the daily light cycle. The rats had ad libitum access to water at all times. The room light intensity during the 12-hour light period ranged from 38–40 foot-candles measured 1.5 m from the floor. The rats as controls from a toxicogenomic study received a dose of continually stirred 0.5% aqueous methylcellulose at 5 ml/kg body weight 4 hour after light on or 4 hour after lights off and were sacrificed at 6, 18, 24, and 48 hours after oral gavage. The necropsies took place within 1 hour of the scheduled period of time. Dosing in the light reversal group was accomplished under a dim red light (<0.2 lux with a wavelength of greater than 650 nm).

The rats were anesthetized with CO₂/O₂, blood samples collected for clinical chemistry by cardiac puncture, the abdominal cavity was opened, and the portal vein severed before necropsy.

Clinical Chemistry

Blood was collected and placed into vacutainers with no additives, allowed to clot for 30 minutes, centrifuged, and serum harvested. The serum samples were analyzed on a centrifugal analyzer (Hitachi 911, Chula Vista, CA). Parameters analyzed included alanine aminotransferase (ALT), aspartate aminotransferase (AST), sorbitol dehydrogenase (SDH), bile salts/acids, alkaline phosphatase, total protein, albumin, blood urea nitrogen, creatinine, and cholesterol. Serum was processed as it was collected in 8 runs conducted over 5 days for the 2 experiments. The serum clinical chemistry results were compared against controls that were included for each run.

Histological Methods

After blood collection, the liver was promptly removed. Half of the left lobe and half of the median lobe was placed in RNALater (Ambion, Austin, TX) for subsequent RNA isolation. A section was taken from the left and median lobes for histology and placed in 10% neutral buffered formalin (NBF). The liver sections were embedded, sectioned, stained with H&E, and examined for pathological changes. The diagnoses were subject to QA and pathology peer review following standard NTP procedure (Boorman et al., 2002).

Hepatic Glutathione Analysis

The method was based on the determination of thiols by high-performance liquid chromatography with fluorescence detection (Toyo’oka et al., 2001). Samples of the left hepatic lobe from individual rats were analyzed for reduced and total glutathione. One mg of liver was placed in a test tube with 1 ml of 5% trichloroacetic acid containing 5mM EDTA, homogenized, vortexed, and centrifuged. Next, 50 μl of homogenate was transferred to a test tube with 950 μl of 5 mM EDTA containing 0.1 M borax (pH 9.3). Then, 100 μL aliquots were placed in duplicate tubes for each liver sample.

Fluorobenzofurazan-4-sulfonic acid ammonium salt was added for reduced glutathione and tri-n-butylphosphine in acetonitrile was added to the other tube for total glutathione. High performance liquid chromatography with fluorescence detection was used for determining total and reduced glutathione and compared against calibration standards that were processed with each assay.

RNA Isolation

The liver tissues collected for RNA isolation were chopped with sharp razors into 0.5-cm cubes or smaller in RNALater within 4 minutes of necropsy. The tissues were stored in RNALater overnight at 4 ±3°C and then stored at −20 ±1°C until RNA isolation (within 60 days). For RNA isolation 130–150 mg of tissue from the left lobe was weighed out, 1 ml of chilled lysis buffer was poured over the tissue that was then minced into 1-mm pieces, 7 ml of lysis buffer was then added and the tissue homogenized with a handheld Omni tissue homogenizer with a disposable plastic 7-mm-diameter Omni generator probe (Omni # 34750, Omni International, Marietta, GA) at no more than half maximum speed for 45 seconds. The tissue lysate was centrifuged for 5 minutes at 4000 × g, the supernatant divided into 2 tubes, 3.5 ml of 50% absolute ethanol added and the tubes were vigorously shaken for 1 minute. RNeasy midi spin columns (Qiagen, Valencia, CA) were used for RNA isolation. The RNA was concentrated using Millipore Microcon centrifugal filter devices (Billerica, MA), analyzed using a spectrophotometer (capable of reading absorbance at 260 nm and 280 nm with UV light), frozen at minutes −70°C and shipped to the NTP repository until transfer to Paradigm Genetics Inc. (Research Triangle Park, NC) for microarray analysis.

Microarray Hybridizations

The study design included two groups of three rats at each time of the day sacrificed in 2 independent experiments (Figure 1). Six rats formed each composite pool with 2 pools at each time of day. From the pool, 3 individual rats (biological replicates) were selected (randomized so as to include 1 or 2 rats from each independent experiment). Then 1 μg of total RNA from either an individual rat or from a pooled sample was amplified and labeled with a fluorescent dye (either Cy3 or Cy5) using Agilent Technologies’ Low RNA Input Linear Amplification labeling kit (Palo Alto, CA) following the manufacture’s protocol. The use of linearly amplified RNA in gene expression profiling studies has been explored and validated previously (Baugh et al., 2001; Pabon et al., 2001). The amount and quality of the fluorescently labeled cRNA was evaluated using a Nanodrop ND-100 spectrophometer (Nanodrop Technologies, Wilmington, DE) and an Agilent Bioanalyzer. Equal amounts of Cy3 and Cy5-lableled cRNA (750 ng) from the individual rat and from the pooled control, respectively, were hybridized to an Agilent Rat Oligo Array and the fluors were reversed for a second (dye reversal) hybridization. Therefore, 48 hybridizations were performed for the 24 individual rats examined.

Data Analysis

Data from dye reversal hybridizations representing the same individuals were combined in the microarray analysis software package Rosetta Resolver version 3.2.2 (Rosetta Biosoftware, Seattle, WA) using a weighted average. From this combined data, genes that were significantly (Rosetta error model, p ≥ 0.01) up- or down-regulated for each individual rat were identified.

When these lists were combined, a total of 8833 features were differentially expressed in at least one rat when compared to its time- and daylight-matched pooled control. Data for these features for all 24 rats were subjected to an agglomerative clustering algorithm in which Ward’s minimum variance (Ward, 1963) was used as the heuristic criteria and Euclidean distance as the similarity metric. The clustering and resulting dendogram produced was performed in Rosetta Resolver.

Evaluations of associations between toxicology parameters (body weights, clinical chemistry results, hepatic glutathione levels, cage rank by weight and experiment) were performed using supervised and unsupervised methods. Principal Components Analysis (PCA) was used to evaluate any global differences between the samples and identify the primary sources of variation. Supervised analysis to find genes associated with toxicological parameters was performed using Significance Analysis of Microarrays (SAM) (Tusher et al., 2001). The SAM false discovery rate was set as less than or equal to 1% for these analyses.

Results

This evaluation was conducted to determine the variation in hepatic gene expression in study control groups. The sacrifices at 6, 18, 24, and 48 hours after dosing with 0.5% aqueous methylcellulose simulate vehicle-control animals in a typical toxicogenomic study and are part of a larger study to evaluate the effect of time on hepatotoxicity. Twenty-four individual rats were analyzed for gene expression from a total of 48 vehicle treated male rats (Figure 1). Clinical chemistry and other toxicology results were generated for the entire set of 48 rats.

Histopathological evaluation of the left hepatic lobes revealed no specific pathological diagnosis for any of the rats. Occasional minimal focal infiltrates of mononuclear cells were found in the liver and were considered within normal range for untreated rats of this age.

The 48 rats had a mean body weight of 286.4 ±18.4 g (mean ±SD) with a range from 252 to 325 g. The mean weight of the subgroup of 24 selected for individual hybridization was 288.5 ±9.9 ranging from 264 g to 314 g. Rats that were 10 g heavier than the smallest rat in the cage were considered alphas and those 10 g less that the alphas were considered omegas (Table 1) for purposed of comparing gene expression by cage rank. The hepatic reduced glutathione levels ranged from 2.1 to 8.0 μM/g of liver. Analysis of variance to determine whether hepatic glutathione levels varied with other biological parameters found no correlation between hepatic glutathione levels and body weight, cage rank, serum activities of ALT, AST, SDH, bile acid concentrations, or between experiments 1 and 2. There was a weak correlation ( p = 0.1) between hepatic glutathione levels and time of day with the highest reduced hepatic glutathione values found in rats sacrificed at noon. Selected clinical chemistry data and other relevant information from the 24 rats selected for individual hybridizations are shown in Table 1.

The two independent studies were compared for differences between the parameters shown in Table 1. As expected there were no differences in body weights, hepatic glutathione levels, or clinical chemistry variables, with the exception of ALT activities, for the 2 experiments. Unexpectedly, serum ALT activities were 62 ±17 IU/L (mean ±SD) in experiment 1 versus 76 ±15 IU/L in experiment 2 ( p = 0.004). The historical control values for ALT at the laboratory are 71 ±17 IU/L (mean ±SD). Two rats with ALT levels of 116 and 120 IU/L fell outside 1 standard deviation of the historical mean. These 2 rats were removed and the data were reanalyzed. The difference between ALT levels in study 1 (59 ±13 IU/L) and study 2 (74 ±12 IU/L) remained highly significant ( p = 0.0003). There was no association between ALT levels and time of day, body weight or cage rank. By chance, 1 of the high ALT rats (ALT of 120 IU/L, rat 140) was selected for individual hybridization and was part of the 6 PM group. The other high ALT rat (116 IU/L, rat 135) was part of one of the RNA pools for rats sacrificed at midnight.

The ALT data and laboratory records were analyzed for possible technical sources of variability. The serum ALT analyses were conducted shortly after necropsy in 8 individual runs over 5 consecutive days from the 2 experiments. There was some overlap between the 2 experiments. On day 3, the 48-hour samples were run for experiment 1 as well as the 6-hour serum samples for experiment 2. Low (expected values 24–32 IU) and high (expected values 92–120 IU) ALT reference controls were run before and after each group of serum samples. The ALT values for the reference controls run with the serum samples were analyzed for systematic bias across the days. There was no evidence of bias in the ALT reference standard values across the days. AST and SDH also used as biomarkers of hepatic damage did not show statistical differences between experiment 1 and experiment 2.

Gene expression profiling was used to determine whether body weight, serum chemistry, or cage rank by weight had any consistent effect on transcript levels. Figure 1 details which individual rats were selected for analysis (number in red) and the six rats that contributed to each control pool. The biological parameters for each rat are shown in Table 1. Genes demonstrating significant differential expression ( p ≥ 0.01) compared to their pooled controls were identified using Rosetta Resolver system (Seattle, WA). When combined, a total of 8,833 genes were up- or down-regulated in at least 1 of the 24 rats examined.

An unsupervised method utilizing all the differentially expressed genes that occurred in at least 1 control animal was to determine if single or multiple biological/clinical chemistry parameters were contributing or indicative of these altered transcript levels. Thus, hierarchical clustering was performed on the 24 rats using the 8,833 differentially expressed genes. Subsequent analysis using additional filtering criteria revealed the same structural dendogram as the one depicted in Figure 2. The same trend was also observed when all 20,000+ features present on the array were used for clustering.

The unsupervised clustering (Figure 2) shows roughly four branches. There was no evidence of clustering by time of day as expected since the pools against which the individual rats were hybridized were matched for time of day. There was no evidence that the light reversal rats clustered differently from the day rats, suggesting that the animal successfully adapted to light reversal as also demonstrated by the serum melatonin levels. Further there were no genes that were differentially expressed when comparing experiment 1 to experiment 2.

The bottom cluster (Figure 2) includes 3 rats (rats numbered 117, 140 and 144) that appear to differ most from all other individuals. This group included the rat with the highest ALT activity, but the other 2 rats had ALT activities of 87 and 84 U/L. There was no correlation with other parameters shown in Table 1 and thus not apparent why these 3 rats appear different. It is important to note that the range of differential expression in Figure 2 is quite narrow. Less that 1% of the genes showed a greater than 2-fold difference. Using Principal Components Analysis (PCA) approximately 60% of the variability lies in the first principal component and 6% in principal component 2; components 3 to 10 are all less than 5%. Rats numbered 117,140, and 144 (most separation in unsupervised clustering) also demonstrate more array variability but the range of variability across all arrays was small (Figure 3).

Hepatic glutathione levels are increased after feeding and vary with the overall nutritional status of the animals. Therefore, a supervised analysis was conducted using rats (N = 12) with lower (<4 μM reduced glutathione per gram liver) and rats (N = 12) with higher (<4 μM) hepatic glutathione levels. Differentially expressed genes were not found with regard to hepatic glutathione levels.

Male Fischer rats kept under standard toxicology study conditions appeared to have fairly consistent hepatic gene expression levels. The toxicological parameters examined did not appear to correlate with the random differences in individual hepatic gene expression levels.

Discussion

This study compared gene expression variation in the livers of 24 individual control inbred rats housed under tightly controlled environmental conditions. RNA was isolated from the livers of 3 individual rats and compared against a pool of RNA from 6 time-matched controls. Standard toxicological parameters including liver histology, clinical chemistry analysis, body weights, and liver glutathione levels were measured in each rat.

Group housed rats may establish a hierarchy that could affect feeding patterns (Stefanski et al., 2001). While these were age-matched controls from a commercial supplier, there was weight variation (up to 14%) between some cage mates. Hepatic glutathione levels increase with eating and decrease with fasting (Jenniskens et al., 2002; Alhamdan and Grimble, 2003; Mariotti et al., 2004). Hepatic glutathione levels did not correlate with cage rank by weight, experiment, necropsy time or body weight using analysis of variance procedures. The slightly increased hepatic glutathione levels at noon were consistent with rats just ending a feeding period. Hepatic glutathione levels reflect the oxidative state of the liver and capacity to respond to injury (Schauer et al., 2004). Hepatic glutathione levels did not correlate with serum ALT levels or with specific genes. It appears that the range of hepatic glutathione levels in our control animals was small and not predictive for differentially expressed genes. Comparison of hepatic gene expression from heavier rats versus lighter rats in a cage failed to identify any consistent differences. This simple separation is probably not adequate to identify stress or feeding patterns. It may be prudent to record the cage and weights of rats used in microarray experiments until more is known about possible cage effects.

The mean serum alanine aminotransferase (ALT) activity in experiment 1 differed from experiment 2. Concern about possible drift in the ALT assays run over 5 consecutive days for the 2 experiments led us to analyze the high and low reference controls included in each run. The ALT values for the reference controls run with each group of serum samples for days 1 and 2 were not significantly different from reference standard values for days 3–5. AST and SDH activity, other markers of hepatic injury did not differ between the 2 experiments. Since the serum ALT half-life is similar to AST, it is unlikely that this indicates an acute injury that is reflected only in the ALT activity. The ALT activities in both experiments were not different from the historical control value (71 ±17 IU/L). The liver appeared similarly normal on histological examination for the 2 experiments and finally there were no differentially expressed genes between the 2 studies, suggesting the difference in ALT activity may be normal variation.

In healthy middle-aged adults, increased ALT or AST activities exceeding 2 times upper level of normal are found in approximately 1% of the serum samples (Downs et al., 2001). Serum ALT activities are used to screen patients and abnormal ALT activities were found in up to 14% of patients in 1 study (Piton et al., 1998). In our study, the rat with the highest ALT level had hepatic transcripts that differed from the other 5 individual control rats in the pool. This suggests that clinical pathology may be useful in screening rats prior to inclusion in a pooled reference sample.

This study suggests that some interindividual variation occurs in inbred rats under study conditions. There is little evidence of a pattern of specific genes. No changes in gene expression were common to even half of the control rats and only 17 gene alterations were found when restricted to 10 out of the 24 controls. The vast majority of the differentially expressed genes show very low-fold change. In this study 95% of the differentially expressed genes exhibited less than 1.5-fold increase or decrease in expression. PCA analysis revealed that most of the variability was in the first principal component, indicating that there is little structure outside of random variability in the data. Finally, the variability for each array was relatively similar across all arrays suggesting that no outlier arrays are present.

The individual variability in this study suggests that the use of paired controls is not ideal since an alteration in transcript level on a 2-color microarray may be due to a change in either the treated animal or the paired control. The use of 6 rats to form a reference pool as was done here dilutes the effect of an outlier animal and reduces variability attributable to the controls. The use of a standard reference may offer advantages but development of a standard reference and demonstration of its utility in toxicogenomic studies remains to be done.

An understanding of biological variability is critical for recognizing subtle changes that may be induced by a toxicant. Microarray studies can be designed to account for variability in control animals. Study design and careful documentation of the study will facilitate large-scale data mining across studies and refine approaches for identification of biologically relevant markers of toxicity.