|Home | About | Journals | Submit | Contact Us | Français|
The use of herbal dietary supplements in the United States is rapidly growing, and it is crucial that the quality and safety of these preparations be ensured. To date, it is still a challenge to determine the mechanisms of toxicity induced by mixtures containing many chemical components, such as herbal dietary supplements. We previously proposed that analyses of the gene expression profiles using microarrays in the livers of rodents treated with herbal dietary supplements is a potentially practical approach for understanding the mechanism of toxicity. In this study, we utilized microarrays to analyze gene expression changes in the livers of male B6C3F1 mice administered Ginkgo biloba leaf extract (GBE) by gavage for 2 years, and to determine pathways and mechanisms associated with GBE treatments. Analysis of 31,802 genes revealed that there were 129, 289, and 2,011 genes significantly changed in the 200, 600, and 2,000mg/kg treatment groups, respectively, when compared with control animals. Drug metabolizing genes were significantly altered in response to GBE treatments. Pathway and network analyses were applied to investigate the gene relationships, functional clustering, and mechanisms involved in GBE exposure. These analyses indicate alteration in the expression of genes coding for drug metabolizing enzymes, the NRF2-mediated oxidative stress response pathway, and the Myc gene-centered network named “cell cycle, cellular movement, and cancer” were found. These results indicate that Ginkgo biloba-related drug metabolizing enzymes may cause herb–drug interactions and contribute to hepatotoxicity. In addition, the outcomes of pathway and network analysis may be used to elucidate the toxic mechanisms of Ginkgo biloba.
Since the U.S. Congress passed the Dietary Supplement Health and Education Act (DSHEA) in 1994, herbal products have been the fastest growing segment of the vitamin, mineral supplements, and herbal products industry in the United States. In 2004, the American Herbal Products Association estimated that there were about 3,000 species of plants in as many as 50,000 different products sold as herbal supplements in the United States (Zurer and Hanson, 2004). St. John's wort, Ginkgo biloba, golden seal, panax ginseng, kava, aloe vera, and mild thistle extract are among the most widely used of these products (Chan and Fu, 2007; Chan et al., 2007; Fu et al., 2008; Guo et al., 2009).
Although it has been reported that a number of herbal dietary supplements cause adverse health effects (Chan and Fu, 2007; Chan et al., 2007; Fu, 2007; Fu et al., 2007, 2008; Gurley et al., 2005, 2007; Hu et al., 2005; Singh, 2005), to date, safety issues concerning potential side effects and toxic contamination of herbal products have not been adequately addressed. Thus, assessment of the safety of herbal plants and herbal dietary supplements is timely and important (FDA, 2001, 2004a, 2004b; Fong, 2002; Fu et al., 2002). Recently, a number of herbal dietary supplements and active ingredients have been nominated by the U.S. Food and Drug Administration (FDA) and the U.S. National Institutes of Health (NIH) to the National Toxicology Program (NTP) for determination of their toxicity and tumorigenicity. Ginkgo biloba, panax ginseng, kava, aloe vera, and green tea are among the herbal dietary supplements currently under investigation by the NTP.
In general, approaches for determining the mechanism by which a pure chemical induces toxicities or tumors have been well established. Nevertheless, it is still a challenge to determine the mechanisms of toxicities or tumor induction elicited by a mixture of many chemical components, such as herbal plants and herbal dietary supplements. For chemical mixtures there is a need for new approaches for elucidating mechanisms. Toward this goal, we previously examined the alterations in gene expression of drug metabolizing enzymes in the livers of Fischer 344 rats administered kava extract by gavage for 14 weeks (Guo et al., 2009). Our results indicate that kava extract can significantly modulate drug-metabolizing enzymes, which can potentially cause herb–drug interactions and may be responsible for different types of liver injuries. The gene expression profile correlated well (Guo et al., 2009) with immunohistochemical data collected using the same liver tissues (Clayton et al., 2007). We observed that kava altered the expression of Cyp1a1 (Cytochrome p450 1a1) and many other Cyp enzymes that metabolize various xenobiotics and drugs, and the changes for many of those genes were first reported by us (Guo et al., 2009). These findings illustrate that, without obtaining the whole spectrum of gene expression changes, some important information may be missed. Our study also suggested that analysis of the gene expression profiles using microarrays in the livers of rodents treated with herbal dietary supplements is a potentially practical approach for understanding the mechanism of toxicity (Guo et al., 2009).
Ginkgo biloba has been one of the most widely sold products in health food stores in the United States, with total sales exceeding $100 million (Chan et al., 2007). Several reports have indicated that Ginkgo biloba inhibits Cyp activity and, when taken in combination with prescription and conventional medications, may produce Cyp-mediated herb–drug interactions (Bressler, 2005; Gurley et al., 2005, 2007; Hu et al., 2005; Matthias et al., 2007; Singh, 2005; Williamson, 2005; Wold et al., 2005; Zou et al., 2002). Gingko biloba leaf extract was nominated by the National Cancer Institute to the NTP to conduct toxicological evaluation, mechanistic studies, and a 2-year chronic carcinogenicity bioassay based on the following reasons: (1) Ginkgo biloba is a well-defined product and itself or its active ingredients have clearly demonstrated biological activities, (2) Ginkgo biloba can be consumed in rather large doses for an extended period of time, and (3) some ingredients in Ginkgo biloba are known mutagens (NTP, 1998). As a continuation of our mechanistic studies of herbal dietary supplements, in the present study, gene expression changes were analyzed with the focus on drug metabolizing genes in the livers of male B6C3F1 mice treated orally and chronically with Ginkgo biloba extract. In addition, changes in pathways and networks, which may indicate toxicity, were also explored.
Ginkgo biloba leaf extract (GBE) (CAS No. 90045-36-6) was received from the Midwest Research Institute (MRI, Kansas City, MO). Infrared spectroscopy and HPLC/UV analyses confirmed the chemical was GBE by comparing with a reference sample received from MRI. The identity, chemical purity, and formulation analysis evaluation data are available at the archives of the National Institute of Environmental Health Sciences. Male B6C3F1 mice (38–44 days old) were obtained from Taconic Laboratory Animals and Services (Germantown, NY) and were used in the NTP toxicity and carcinogenesis studies of GBE in mice. The animals were quarantined for 14 days before separating at random into four groups. GBE in corn oil was administered to each group by gavage at 0, 200, 600, or 2,000mg/kg, 5 days a week for 104 weeks. The animals were monitored twice daily during the course of the experiment, and body weights and clinical observations were recorded weekly for the first 13 weeks and every 4 weeks thereafter, and at the study termination. Animal handling and husbandry were conducted in accordance with guidelines of the National Institutes of Health. During necropsy at terminal sacrifice (104 weeks); a piece of the left lateral lobe of the liver was collected from 10 animals from each group, cut into three pieces and frozen in liquid nitrogen until use.
Total RNA was isolated from liver tissues of control and GBE-treated mice using an RNeasy kit (Qiagen, Valencia, CA). The RNA was quantified by optical density reading using a NanoDrop ND-1000 spectrophotometer (NanoDrop products, Wilmington, DE). The purity and quality of extracted RNA were evaluated using the RNA 6000 LabChip and Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA). RNA samples with RNA integrity numbers (RINs) greater than 9.0 were used for target labeling, microarray experiments, and TaqMan assays.
Five micrograms of RNA from each sample was used for Cy5-labeled aRNA preparation using the Ambion MessageAmp aRNA Amplification Kit (Applied Biosystems Inc, Foster City, CA) with minor modifications to adopt the amino-allyl UTP indirect labeling protocol. In brief, the 4μL T7 UTP solution (75mM) was replaced with 2μL of T7 UTP solution (75mM) and 3μL of 50mM 5-(3-aminoallyl)-UTP in the process of in vitro transcription to synthesize aRNA. The reaction conditions and the purification of unlabeled aRNA remained the same as described in the user's manual. For Cy5 labeling of aRNA, 40μg of unlabeled aRNA was used in each labeling reaction to generate sufficient Cy5-labeled cRNA for the microarray hybridization in triplicate. The Cy5 dye (GE Healthcare, Piscataway, NJ) was suspended in 32μL DMSO before use. Eight microliters of Cy5 dye solution and 2.5μL of 10μM NaOH were added to the 40μg aRNA and the volume was adjusted to 25μL by addition of ddH2O. The mixture was incubated at 25°C in the dark for 2h. The labeled cRNA was then purified using the Ambion NucAway Spin Column (ABI). The amount and labeling efficiency of the purified Cy5-labeled cRNA were quantified by optical density reading using a NanoDrop ND-1000. The minimal requirement of labeling efficiency was 10 Cy5 dye molecules per 1,000 nucleotides.
Microarray experiments were performed using Phalanx Mouse OneArray Version 1.1 (MOA 1.1; Phalanx Biotech Group, Inc., HsinChu, Taiwan). Each microarray contains 31,802 oligonucleotide probes that include 29,922 mouse gene probes for transcription expression profiling and 1,80 experimental control probes. The MOA 1.1 content is based on an abridged version of the Mouse Exonic Evidence Based Oligonucleotide (MEEBO). MEEBO consists of a set of 70-mer probes made available to the public specifically for DNA microarray. For microarray hybridization experiments, 10-μg Cy5-labeled cRNA derived from each RNA sample was applied to each of three MOA 1.1 microarrays. Prior to the microarray hybridization, the Cy5-labeled cRNAs were fragmented using the reagents and protocol provided in Ambion RNA Fragmentation Reagents kit (Ambion Inc., Austin, TX). Each 10-μg fragmented Cy5-labeled cRNA was suspended in OneArray hybridization buffer at a final volume of 180μL for microarray hybridization. After microarray hybridization, the arrays were scanned using a fluorescence scanner (GenePix 4000B; Molecular Devices, Sunnyvale, CA) and the fluorescent intensities were extracted from the generated images by following the instructions and the conditions described in the MOA 1.1 User's Manual.
Each Cy5-labeled cRNA derived from a liver RNA sample was hybridized to MOA 1.1 microarrays in triplicate. The raw microarray intensity data was imported to into ArrayTrack (www.fda.gov/nctr/science/centers/toxicoinformatics/ArrayTrack/). The intensities of the three technical replicates for each probe were averaged and then normalized by quantile normalization. Differentially expressed genes (DEGs) were selected based on t-test and fold change cutoff criteria. Principal component analysis (PCA) and hierarchical cluster analysis (HCA) were conducted within ArrayTrack. Additional calculations were performed within JMP 7.0 (SAS Institute, Cary, NC). The pathways, networks,and functional analyses were generated through the use of ingenuity pathways analysis (IPA;Ingenuity Systems, www.ingenuity.com). Canonical pathways analysis identified the pathways from the IPA library of canonical pathways that were most significant to the data set. DEGs that associated with a canonical pathway in the IPA Base were considered for the analysis. The significance of the association between the data set and the canonical pathway was measured in two ways: (1) a ratio of the number of genes from the data set that map to the pathway divided by the total number of genes that map to the canonical pathway is displayed, and (2) Fischer's exact test was used to calculate a p-value determining the probability that the association between the genes in the dataset and the canonical pathway is explained by chance alone. DEGs containing gene identifiers and corresponding expression values were uploaded into the application to generate networks. Each gene identifier was mapped to its corresponding gene object in the ingenuity pathways knowledge base. These genes, called focus genes, were overlaid onto a global molecular network developed from information contained in the ingenuity pathways knowledge base. Networks of these focus genes were then algorithmically generated based on their connectivity.
DEGs identified with microarray analysis were selectively confirmed by TaqMan assays (Applied Biosystems, Foster City, CA). Eleven TaqMan probes were used in these assays, including Cyp1a1 (Mm00487218_m1), Cyp1a2 (Mm00487224_m1), Cyp2b10 (Mm00456591_m1), Cp2b13 (Mm00771172_g1), Cyp2c55 (Mm00472168_m1), Cyp2d13 (Mm00775259_gl), Ces2 (Mm00524035_m1), Fmo3 (Mm00514964_m1), Gsta2 (Mm00833353_mH), Actb (Mm00607939_s1), and Polr2a (Mm 00839493_m1). Actb and Polr2a were used as endogenous controls. cDNA was prepared using a High-Capacity cDNA Archive Kit (Applied Biosystems), that is, 2μg RNA was reverse-transcribed in a final volume of 20μL with random primers at 25°C for 10min followed by 37°C for 120min according to the manufacturer's instructions. Each TaqMan assay was run in triplicate for each RNA sample. Total cDNA (25ng) in a 25μL final volume was used for each assay. Assays were run with Universal Master Mix (2×) without AmpErase UNG on an Applied Biosystems 7900 HT Real-Time PCR System using universal cycling conditions (10min at 95°C; 15s at 95°C, 1min 60°C, 40 cycles).
Two endogenous control genes, β-actin (Actb) and RNA polymerase II A (Polr2a), were used for normalization. Each replicate cycle threshold (CT) was normalized to the average CT of the two endogenous controls on a per sample basis. The comparative CT method (Livak and Schmittgen, 2001) was used to calculate relative quantification of gene expression. The following formula was used to calculate the relative amount of the transcripts in the GBE-treated sample (treated) and the vehicle-treated sample (control), and both were normalized to the endogenous controls: ΔΔCT=ΔCT (treated)−ΔCT (control), where ΔCT is the difference in CT between the target gene and endogenous controls by subtracting the average CT of controls. The fold change for each treated sample relative to the control sample=2−ΔΔCT. A list of DEGs was identified using a two-tailed t-test. The criteria were P value less than 0.05 and a mean difference equal to or greater than two-fold.
During the course of the study no clinical findings attributed to GBE administration were observed. After 18–19 months of GBE treatment, the body weights of the 600 and 2,000mg/kg mice began to decrease gradually. At terminal sacrifice the body weights of these two dose groups of mice were significantly lower (~17 and ~25%, respectively) compared with controls.
In this study, the global changes of gene expression in mouse liver treated with GBE were determined by high-density oligonucleotide microarray analysis (Mouse OneArray, which contains 31,802 probes covering 29,922 verified mouse genes). Total RNA was isolated from liver tissues of control mice and mice treated with 200, 600, and 2,000mg/kg GBE. Each group consisted of four animals; therefore, 16 RNA samples were isolated for microarray analysis. Three technical replicate arrays were performed for each RNA sample, and thus a total of 48 microarrays were used in this study.
The overall quality of data from the 48 microarrays was first assessed by HCA. The log2 intensity of the entire gene set including 31,802 probes was scaled by Z-score transformation, and then these values were hierarchically clustered using a distance metric of 1−R, where R is the Pearson correlation coefficient between two samples, and the average linkage (Fig. 1). As expected, for most samples, the three technical replicates for the same RNA sample were tightly grouped, indicating close similarity among them. In addition, two apparent clusters were revealed. The first cluster (in black and red) consists of samples from the control, and low- (200mg/kg) and middle-dose (600mg/kg) GBE treatments; the second cluster (in blue) includes samples from high-dose (2,000mg/kg) GBE treatments. Samples in the high-dose group were clearly separated from those in the control, low, and middle groups (Fig. 1). In the first big cluster, the three subgroups (control, low, and middle dose) were also well separated with few exceptions.
The reproducibility of the three technical replicates was further assessed by calculating the Pearson correlation coefficient of pairwise log2 intensity. R values of each pairwise comparison were obtained; the median R values of the three technical replicates are shown in Table 1. The median R values across 16 RNA samples range from 0.984–0.997. One array from sample C38 (C38_3) and one from C24 (C24_2) showed dissimilarity in comparison with its other two technical replicates with relatively low correlation coefficients (0.979 and 0.965, respectively). These two hybridizations were considered as outliers and were removed from further calculation and data analysis.
It is worth noting that the intensities of all data points (31,802 probes) were included for the comparison. No specific cutoff was applied for evaluation of microarray quality. For visualization purposes, as an example, the raw log2 intensity data of all 31,802 probes from the three technical replicates (C20_0_1, C20_0_2, and C20_0_3) of animal #C20 were plotted against each other (Fig. 2). As shown in the scatter plot, for most spots, the intensity values from the replicate microarrays accumulated along the diagonal axis with a correlation >0.996, indicating that the data from the technical microarrays were highly reproducible.
The intensity values of each probe from the three technical replicates were averaged and log2 transformed. The overall quality of the microarray data from the 16 RNA samples was further assessed by calculating the Pearson correlation coefficients of the averaged log2 intensity data for all pairwise sample comparisons within each control or treatment group (Table 2). The median correlation coefficient was 0.973 across all 16 arrays with a range from 0.926 to 0.989. The median correlation for the control group was 0.987, with a range from 0.983 to 0.988, whereas the median correlation for the 2,000mg/kg treatment group was 0.970, with a range from 0.961 to 0.977. These results demonstrated the high degree of gene expression similarity for samples within the same group. On the other hand, the median correlation between the control group and the 2,000mg/kg GBE treated group was 0.957, with a range from 0.933 to 0.964, demonstrating that there were significant differences in the gene expression profiles between the control and high dose GBE treated groups.
The log2 transformed intensity of any two gene expression profiles was plotted and compared, as shown in Figure 3. When intensities of the entire probe set obtained from a control sample (C20) were plotted against another control sample (C36), most of data points gathered along the diagonal axis of the scatter plot with a correlation coefficient value of 0.987, demonstrating the good repeatability of the two biological replicates (Fig. 3A). However, in comparison between the control sample (C20) and the 2,000mg/kg GBE-treated sample (C32), many data points were scattered, indicating a great number of genes were altered in response to the treatment of high-dose GBE (Fig. 3B).
The intensities of the gene expression data were also analyzed by PCA. PCA uses analysis of the principal sources of variance in data and displays this information graphically either in a two-dimensional or three-dimensional space, which is most commonly used to classify gene expression profiles (Wang and Gehan, 2005). Figure 4 displays a PCA 3D view using the first three principal components for gene expression profiles from the samples. The PCA result showed a clear separation between controls and 2,000mg/kg-treated samples, suggesting that there was a treatment effect on liver gene expression which is consistent with the results from the HCA (Fig. 1). Less separation was observed between controls, low-dose (i.e., 200mg/kg), and middle-dose (i.e., 600mg/kg) groups based on PCA.
The average intensities for each probe from three technical replicates of each sample were normalized by the quantile normalization method, which assumes that a quantile plot of two data vectors with the same distribution will have a straight diagonal line. The normalized data were used for selection of DEGs. We had reported that a straightforward approach of fold change selection plus a nonstringent P cutoff is useful in identifying reproducible gene lists (Guo et al., 2006; Shi et al., 2006). We used a twofold change in gene expression compared to the controls and a P-value less than 0.01 for the difference as minimum requirements for the selection of DEGs. Based on these two criteria a total of 129, 289, and 2,011genes, respectively, in 200, 600, and 2,000mg/kg treatment groups satisfied the requirements and were identified as DEGs. Out of 2,011genes identified in the 2,000mg/kg treatment group, 817 were upregulated and 1,194 genes were downregulated in response to GBE exposure. Significance analysis of microarray (SAM) (Tusher et al., 2001), the statistical method that calculates difference in gene expression based on permutation analysis of expression data and calculates a false discovery rate, was also applied to generate DEGs. With an FDR of 0.01 and fold change greater than 2, a total of 2,211genes were identified in the 2,000mg/kg treatment group compared to control group. Out of 2,011genes that were indentified as DEGs using the t-test, 1,952genes (97%) are overlapped with those produced by the SAM, indicating that similar gene lists were produced with two gene selection methods.
Because metabolic activation of chemicals is very important for liver toxicity, and alteration of drug metabolizing genes is suspected to contribute to the toxicity of Gingko biloba, we first focused our investigation on the gene expression changes of drug metabolizing enzymes and found that the expression of a large number of drug-metabolizing genes was altered. Table 3 lists the 68 drug-metabolizing genes whose expression was significantly altered by the 2,000mg/kg GBE treatment. Among the 68 altered genes, 36genes were upregulated and 32genes were downregulated. As tabulated in Table 3, among the 68 drug-metabolizing enzyme associated genes, 33genes are Phase I metabolizing enzymes; 18genes are Phase II metabolizing enzymes, and 17genes are transporters (Phase III). As summarized in Table 4, 21% (33 out of 154) of Phase I genes, 21% (18 out of 87) of Phase II genes, and 24% (17 out of 72) of Phase III genes were altered. Although a relatively small portion of genes (6.7%) were altered at the genome level (2011 DEGs out of 29,922 probes), the percentage of changed drug metabolizing genes (22%, Table 4) was high.
TaqMan assays were used to verify the results for a selected group of genes whose expression changes were seen using microarrays. Nine drug metabolizing genes were selected for the TaqMan validation. The results are shown in Figure 5A. Based on triplicate measurements for each RNA sample, the genes Cyp1a1, Cyp1a2, Cyp2b10, Cyp2c55, Ces2, and Gsta2 were overexpressed, and the changes in Cyp1a1, Ces2, and Gsta2 were dose dependent. The expression of Cyp2b13, Cyp2d13, and Fmo3 genes were decreased, notably for the highest dose treatment (Fig. 5B). Over 100-fold decrease in gene expression was observed for the three downregulated genes in the 2,000mg/kg treatment group. As expected, TaqMan data were in good agreement with the microarray data (Fig. 5 and Table 3).
IPA (version 7.0) was used to determine the most relevant biological functions, pathways, and networks of the genes altered by 2,000mg/kg GBE treatment. Out of 2,011 DEGs, 1,246 were identified by IPA and were overlaid onto a global molecular network developed from information contained in the database. Not to our surprise, the top canonical pathway was “metabolism of xenobiotics by cytochrome P450” with 36genes involved and a P-value of 1.26E-11 followed by “fatty acid metabolism” (Table 5 and Table 6). Interestingly, NRF2 (nuclear factor erythroid-related factor 2)-mediated oxidative stress response was identified with 23genes involved and a P-value of 5.37E-03. In the pathway regulated by NRF2, most of the genes were upregulated, whereas two genes were downregulated (Fig. 6). The upregulated genes included those coding for detoxifying enzymes (Nqo1, Gstms, Gstas, Gstp1, and Ugt1a2), glutathione homeostasis (Gsr), and transporter (Abcb1). Maf2, a transcriptional factor forms heterodimers with Nrf2, was also upregulated. The two downregulated genes were P13k, which plays a role in Nrf2 phosphorylation, and Sod3, which is a family member of antioxidants. Nqo1 showed a prominent change in expression, exhibiting 7.5- and 11-fold induction for 600 and 2,000mg/kg treatments, respectively.
The significantly changed genes after GBE administration were mapped into 87 networks. Each network was associated with specific genes and involved different functions. Twenty-three networks had 20 or more focus molecules. The first network, which contains the largest number of genes (34genes), with functions related to “cell cycle, cellular movement and cancer,” is illustrated in Figure 7. The central gene in this network is Myc, the pro-oncogene, which was upregulated by about sixfold. Twenty-one genes were found to be associated directly with Myc and these genes function in cell cycle, cell growth/proliferation, and cell death processes. Out of these 21 Myc-associated genes, 15 were upregulated and 6 were downregulated. Notably, Cidec (cell death-inducing DFFA-like effectors c), which promotes apoptosis, exhibited a 47-fold induction. Known molecules of DNA replication checkpoint (Cdca8 and Cdc45l) were also upregulated. Gene Gas1 (growth arrest-specific 1), the tumor suppressor gene that blocks entry to S phase, and prevents cycling of cells, was remarkably downregulated by 20-fold.
Currently there is no established methodology for determining the mechanisms of toxicity (particularly tumorigenicity) induced by a mixture containing many chemical components such as GBE and other herbal plant extracts. Previously, we have proposed that DNA microarray analysis should be a highly practical initial approach for revealing the whole spectrum of gene expression alterations by a chemical mixture (Guo et al., 2009). Microarray technology is a useful tool to rapidly detect the induction/inhibition of drug-metabolizing enzymes after toxicant treatment, and it has greatly contributed to the understanding of drug metabolizing enzyme function and expression (Blomme et al., 2009; Rezen et al., 2007). Although levels of gene expression do not fully represent the levels of enzyme activities, investigations at the gene expression level have revealed that there is a high degree of correlation for Phase I enzymes between the fold inductions of the enzymatic activity and mRNA expression in liver samples (Iyer and Sinz, 1999; Roymans et al., 2004).
In this study, we analyzed gene expression changes in the livers of B6C3F1 mice treated chronically with 200, 600, and 2,000mg/kg of GBE for 2 years using a genome-wide microarray approach. A total of 31,802 probes covering 29,922 verified mouse genes were analyzed, and three technical replicate arrays were performed for each RNA samples. The results of HCA (Fig. 1), Pearson's correlation coefficient of pairwise log2 intensity correlation (Table 1), and scatter plotting of the overall gene expression profiles (Fig. 2) demonstrated the high quality and reproducibility of DNA microarray data obtained from the three technical replicates.
When the gene expression results of GBE-treated rats were compared to the controls, a total of 2,011genes in the 2,000mg/kg treatment group were differentially expressed, of which there were 68 drug-metabolizing enzyme-associated genes, with 33, 18, and 17genes belonging to Phase I, Phase II, and Phase III, respectively (Table 3). It is worth noting that 23 out of the 33 expressed Phase I metabolism associated genes were the Cyps, in which 2genes belong to the Cyp1 superfamily (Cyp1a1 and 1a2), 15genes belong to the Cyp2 superfamily, and 2genes belong to the Cyp3 superfamily (Table 3). The Cyp1, Cyp2, and Cyp3 superfamilies are the most important metabolizing enzymes in the metabolism of drugs and the metabolic activation of toxic and carcinogenic xenobiotics (Gonzalez and Gelboin, 1994; Gonzalez and Yu, 2006).
One of the most reported adverse effects caused by herbal products is hepatotoxicity. Although the mechanism of hepatotoxicity initiated by herbal products is not clear, modulation of drug-metabolizing enzymes, leading to herb–drug interactions, is probably one of the major causes (Clayton et al., 2007; Clouatre, 2004; Guo et al., 2009). The alteration of a number of drug-metabolizing enzymes by GBE treatment (Table 3) could potentially lead to herb–drug interactions as well as changes in the metabolic activation of carcinogenic chemicals when concomitantly present in the body. Our results are consistent with the report that of Ginkgo biloba modulated Cyp 450 enzyme activity, and consequently may produce Cyp-mediated herb–drug interactions (Dubey et al., 2004; Gurley et al., 2004, 2005; Zou et al., 2002).
It is well established that the Cyp1 superfamily is important in the metabolism of xenobiotics, and the members of Cyp2 and Cyp3 superfamilies catalyze metabolism of drugs and other substances. Cyp1a1, which generally exists in a very low amount in the adult liver (Martignoni et al., 2006), was dramatically induced by GBE treatment, as detected by both the microarray analysis and TaqMan assay (Table 3 and Fig. 3A). Cyp1a1 and Cyp1a2 are regulated by AHR, the aryl hydrocarbon receptor, which plays an important role in cell cycle regulation, apoptosis, and development (Xu et al., 2005). It is evident that some toxins such as tumorigenic polycyclic aromatic hydrocarbons (PAHs) (Martignoni et al., 2006) and carcinogen 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) are able to induce Cyp1a1 through AHR, and the upregulation of Cyp1a1 is considered the indicator of TCDD-induced carcinogenesis (Brunnberg et al., 2006). Thus, induction of Cyp1a1 may contribute to the toxic mechanism of Ginkgo biloba, and/or may enhance metabolic activation of the hazardous tumorigenic xenobiotics and jeopardize health.
Besides Cyp, the flavin-containing monooxygenase (Fmo) is another major group of Phase I drug-metabolizing enzymes catalyzing the oxidative biotransformation of drugs. Fmo enzymes break down compounds that contain nitrogen, sulfur, or phosphorus. Fmo3, a major Fmo isoform, plays an important role in metabolizing drugs such as the anticancer drug tamoxifen, the pain medication codeine, the antifungal drug ketoconazole, and antipsychotic drugs clozapine and olanapine (Ciraulo et al., 2005; Klick and Hines, 2007; Parte and Kupfer, 2005). In our study, GBE inhibited the expression of Fmo3 remarkably with the high dose treatment (77-fold reduction detected by microarray and 166-fold reduction validated by TaqMan assay). These results suggest that downregulated expression of the Fmo3 enzyme by GBE potentially inhibits metabolism of certain drugs resulting in drug–drug or drug–food interactions, thus reducing the therapeutic effects or enhancing the toxicity induced by their metabolites.
It is worth noting that although the microarray studies on the changes of drug-metabolizing enzyme expression in kava-treated rats (Guo et al., 2009) and kava-treated mice (Guo et al., submitted for publication) used a 90-day treatment of kava, the present GBE study was performed over a 2-year period of treatment. The difference in time of herbal treatment suggests that the effect of such a dietary supplement on drug metabolizing enzymes is a persistent phenomenon (up to 2 years). Our previous study with kava and the present study with Ginkgo biloba provide substantial evidence that the two herbal dietary supplements cause remarkable changes in drug metabolizing enzymes, particularly the Cyp isozymes, raising the possibility that these herbal supplements might profoundly affect the pharmacokinetics of many coadministrated drugs or other food supplements, thus leading to hepatotoxicity.
In addition to the investigation of alteration for many drug metabolizing genes, in this study we also explored pathways and networks in response to GBE treatment and found that Ginkgo biloba exposure resulted in the significant stimulation of the Nrf2-mediated oxidative pathway (Fig. 6) and the alteration of the Myc-centered network (Fig. 7).
A previous study found that long-term exposure to kava could result in the perturbation of the Nrf2-mediated pathway and eventually generate reactive oxygen species (ROS). In this study, Ginkgo biloba treatments also stimulated the Nrf2-mediated oxidative stress response pathway (p=0.005). Keap1–Nrf2–ARE signaling plays a critical role in protecting cells from endogenous and exogenous stresses, and is involved in antioxidative response, detoxification of xenobiotics, and proteome maintenance (Kensler et al., 2007). Under normal physiological conditions the transcription factor Nrf2 localizes in the cytoplasm and interacts with Keap1. Upon oxidative stress, Nrf2 is released from Keap1 and translocates to the nucleus and subsequently activates its various downstream target genes (Kensler et al., 2007). The target genes show a wide spectrum of functions, such as inactivating oxidants, increasing the levels of glutathione, and enhancing toxin export via transporters to enhance cell survival. As illustrated in Figure 6, in response to GBE treatment, Nrf2 regulated genes including various Gst genes which are involved in GSH-conjugate formation were altered. The Gst genes (Gsta2, Gsta4, Gstms, and Gstp1) and Ugt1a2 were all upregulated following GBE treatment. In turn, the resulting enhancement of GST enzymes is used for neutralizing the electrophiles, the process considered generally as an important detoxification mechanism. NADPH-dependent enzyme NQO1 was also upregulated, which has protective roles toward detoxification of xenobiotic carbonyls and quinones. GBE treatment also elevated Gsr (glutathione reductase) gene expression, and GSR regulates cellular GSH homeostasis by catalyzing the reduction of GSSG to GSH using NADPH as a reducing cofactor (Harvey et al., 2009) In addition, GBE treatment reduced the activity of Sod3, which is one of the superoxide dismutases (SODs) that are the most important line of antioxidant enzyme defense systems against ROS and particularly superoxide anion radicals (Zelko et al., 2002). Studies using Nrf2 knock-out mice showed that they are more susceptible to acetaminophen-induced hepatocellular injury (Chan et al., 2001) and benzo[a]pyrene-induced tumor formation with higher levels of DNA adducts (Ramos-Gomez et al., 2001). This susceptibility is partly due to a reduced level in the expression of detoxification enzymes (Aleksunes and Manautou, 2007; Kensler et al., 2007). Activation of detoxification enzymes plays a pivotal role in protecting cells from oxidative insult when cells encounter a toxin challenge.
Based on the network analysis, 34genes were incorporated into the first network in response to GBE treatment (Fig. 7). The key molecule in this most significantly changed network was Myc pro-oncogene with 5.7-fold increased expression. The protein encoded by this gene is a multifunctional nuclear phosphoprotein that plays a role in cell cycle progression, apoptosis, and cellular transformation. It functions as a transcription factor that regulates transcription of specific target genes. It has been established that stimulation of Myc pathways can cause liver damage, leading to severe lesions up to the tumorigenic level. Overexpression of Myc has been frequently observed in tumors and is the early and critical event in the pathogenesis of hepatocellular carcinoma (Beer et al., 2008; Dang, 1999; Pelengaris et al., 2002). A study using a transgenic model indicated that overexpression of the c-Myc gene sufficiently induced hepatic proliferation and tumorigenesis (Beer et al., 2008; Sandgren et al., 1989). Many studies have shown that Myc gene overexpression during hepatocarcinogenesis induced by carcinogens including choline-devoid diet, arsenite, 5-diethoxycarbonyl-1,4-dihydrocollidine (DDC), and carbon tetrachloride (CCl4) (Beer et al., 2008; Chandar et al., 1989; Chen et al., 2001). In addition, DDC and carbon tetrachloride cooperate with Myc to promote hepatocyte proliferation and rapidly uncover the onset of liver cancer (Beer et al., 2008). Using a toxicogenomic approach to study the molecular mechanisms of nongenotoxic carcinogenicity revealed that major pathways linked to cancer were interconnected via Myc (Nie et al., 2006), and the high percentage of gene expression signature linked to Myc suggested the important role of the Myc gene. Therefore, the alteration of Myc-mediated pathways may serve as potential predictors in carcinogenicity. Among 34genes incorporated into the Myc-centered network (Fig. 7), over 60% of the genes (21genes) were reported to be directly associated with Myc and most of these genes function in cell cycle, cell growth/proliferation, and cell death processes. The Myc associated gene, Gas1 (growth arrest-specific 1) was found to be 20-fold downregulated. Gas1, the tumor suppressor gene, which plays role in growth suppression by blocking entry to S phase, is directly inhibited at the transcriptional level by Myc (Lee et al., 1997), and the transcriptional repression of growth arrest genes is one step in the promotion of cell growth. Downregulation was also one of the distinct events in the process of EGF (epidermal growth factor)-induced hepatocarcinogenesis (Borlak et al., 2005). The loss of growth control through Gas1 may be a necessary event in the multi-step neoplastic transformation.
In summary, we analyzed gene expression profiles in the liver of mice treated with Ginkgo biloba, a widely sold product in health food stores in the United States. Using the systematic approaches of HCA, PCA, pathway analysis, and network analysis, we identified 129, 289, and 2,011genes that were differentially expressed in the livers of 300, 600, and 2,000mg/kg GBE-treated mice, respectively, and observed that a significant portion of the changed genes were associated with drug metabolism, Nrf2-mediated oxidative stress, and the Myc gene-centered network. Based on our results, we speculate that long-term exposure of Ginkgo biloba may cause hepatic damage, including severe liver lesions and even liver tumors. Further studies are warranted.
We thank Drs. Tucker A. Patterson and James C. Fuscoe for their critical review of this manuscript. This article is not an official guidance or policy statement of the National Toxicology Program (NTP) or U.S. Food and Drug Administration (FDA). No official support or endorsement by the FDA or NTP is intended or should be inferred. This research was supported in part by the Intramural Research Program of the NIH, National Institute of Environmental Health Sciences.
The authors declare that there is no conflict of interest.