In Vitro Global Gene Expression Analyses Support the Ethnopharmacological Use of Achyranthes aspera

Achyranthes aspera (family Amaranthaceae) is known for its anticancer properties. We have systematically validated the in vitro and in vivo anticancer properties of this plant. However, we do not know its mode of action. Global gene expression analyses may help decipher its mode of action. In the absence of identified active molecules, we believe this is the best approach to discover the mode of action of natural products with known medicinal properties. We exposed human pancreatic cancer cell line MiaPaCa-2 (CRL-1420) to 34 μg/mL of LE for 24, 48, and 72 hours. Gene expression analyses were performed using whole human genome microarrays (Agilent Technologies, USA). In our analyses, 82 (54/28) genes passed the quality control parameter, set at FDR ≤ 0.01 and FC of ≥±2. LE predominantly affected pathways of immune response, metabolism, development, gene expression regulation, cell adhesion, cystic fibrosis transmembrane conductance regulation (CFTR), and chemotaxis (MetaCore tool (Thomson Reuters, NY)). Disease biomarker enrichment analysis identified LE regulated genes involved in Vasculitis—inflammation of blood vessels. Arthritis and pancreatitis are two of many etiologies for vasculitis. The outcome of disease network analysis supports the medicinal use of A. aspera, viz, to stop bleeding, as a cure for pancreatic cancer, as an antiarthritic medication, and so forth.


Introduction
Large sections of the world population depend on medicinal plants (natural products) for their health care needs. The indigenous population use crude preparations of natural products in the form of decoction, paste, or as food supplement. The prevailing belief is that crude mixtures contain many compounds. Several of these compounds might be inactive by themselves but are essential for the effectiveness of the "active molecule" in the crude preparations. In other words, synergism negates the additive effect of individual molecules [1]. Global gene expression changes effected by nonstandardized (crude) preparations as used by communities at large may not be directly attributed solely to the "active molecule. " Understandably, this approach will not facilitate deciphering the mechanism of action of the active molecules in a crude preparation. However, global gene expression analyses will give a snapshot of the overall effect of crude extract on multiple cellular pathways [2]. Such approaches may explain the possible mode of action of herbal preparations as used in indigenous medicine.
Modern scientific standards emphasize thorough understanding of the molecular mechanism of action of single agent drugs. However, this approach is limited to single molecules, which have been chemically characterized. Even these drugs modulate the activities of disparate metabolic/signaling 2 Evidence-Based Complementary and Alternative Medicine pathways [3]. These results suggest that a drug with specific target can have wide ranging physiological effects.
Analytical methods like microarray, pathway focused PCR arrays, RNA seq, and cGH array give a snapshot of large scale changes at genome level. These methodologies generate large data sets. The development of data mining software has enabled discovery of molecular markers from the mass of data possible. Thus, the ability to analyze gene expression changes at genome level had revolutionized the molecular investigations that until now were difficult. For example, understanding global gene expression changes that accompany certain diseases like cancer has enabled the development of personalized medicine. Data analysis of microarray results also facilitated deciphering molecular pathways that are affected by both existing as well as newly discovered drugs [2].
Achyranthes aspera (Family Amaranthaceae) is a medicinal plant in Ayurvedic pharmacopoeia. It is an annual shrub found in India, Southeast Asia, America, and Sub-Saharan Africa. It is classified as a weed [5]. Traditionally, the plant extract is used as a cure for arthritis, as blood clotting agent, to induce labor (fetal expulsion), and as an abortifacient [6]. Earlier, we demonstrated that methanol-extract of A. aspera leaves (LE) has antiproliferative activity against a variety of human cancer cell lines cultured in vitro. Our results also showed pancreatic cancer cells to be the most sensitive to LE. These exciting outcomes indicate the presence of active anticancer molecule(s) in A. aspera. Unpublished, anecdotal observations suggest that A. aspera benefited pancreatic cancer patients. Animal studies carried out with A. aspera extract confirm LE to be nontoxic to mice. Most importantly, these studies have shown that LE inhibited human tumor growth in mice [7]. In preliminary gene expression studies, we demonstrated LE differentially modulated the expression of key genes involved in vasculogenesis (VEGFs), metastasis (MMPs, TIMPS), apoptosis (Caspase), cell proliferation (Akt-1), and so forth [7,8]. These preliminary results support the idea that LE regulates the function of specific signaling molecules. Except for the above two reports, no additional information is available regarding the global molecular changes caused by LE in human cancer cells. Clearly, simultaneous analysis of global gene expression changes effected by LE would shine light on different pathways affected. Such analyses would help select key genes that can be used as response markers and predict its effect on different illnesses. Based on this premise, we used Agilent Whole Human Genome microarrays to examine the global gene expression changes effected by LE on cultured human pancreatic cancer cell line MIA PaCa-2 (CRL-1420) that were treated with LE for 24, 48, and 72 hours. In this paper, we summarize the effect of methanol extract (LE) of dried powdered leaves of A. aspera on global gene expression pattern in human pancreatic cancer cell line MIA PaCa-2. Using MetaCore software, we narrowed down the pathways affected in MIA PaCa-2 cells treated with LE for three different time points.

Materials and Methods
We have recently published the details about the collection, authentication of Achyranthes aspera leaves, and the preparation of whole leaf extract [7,8]. Briefly, A. aspera leaf powder was serially extracted in 100% hexane followed by 100% acetone over night at room temperature. The acetone extract was concentrated in a rotary evaporator. The residue was dissolved in 100% methanol and called leaf extract (LE).

Cell Lines and Culture
Conditions. The human pancreatic cancer cell line MIA PaCa-2 (CRL-1420) was purchased from ATCC (USA) and maintained in our laboratory in RPMI 1640 medium supplemented with 10% FBS (Invitrogen, USA) at 37 ∘ C in a humidified 5% CO 2 atmosphere. It was seeded in a six-well plate at an initial density of 5 × 10 5 in 2 mL medium. Overnight cultures were exposed to 34 g/mL of LE for 24, 48, and 72 hours, respectively. Two percent methanol exposed cells served as a control. Control and experimental groups were set up in duplicate.

Isolation and Quality
Control of Total RNA. Total RNA was prepared from control and LE treated MiaPaCa-2 cells using TriReagent (Sigma, USA). At the end of the treatment, cells were directly lysed in TriReagent and total RNA isolated. From this point, the rest of the gene expression related procedures were performed at the microarray and gene expression Core facility of the University of Miami. The quality of the extracted RNA was evaluated using BioAnalyzer 2100 (Agilent Technologies, CA). RNA samples with RNA integrity number (RIN) ≥9 have been used for the microarray experiments [9].

cDNA Synthesis and RNA Labeling for Microarray Analyses.
Target RNA was amplified using Amino Allyl Message Amp kit (Ambion, USA). Amplified RNA from the control group was labeled with Cy-3 and the test group with Cy-5 dye (Amersham, USA). Duplicate microarray hybridization was set up using biological replicates that were labeled by dye-swap; that is, control RNA was labeled with Cy-5 and test RNA with Cy-3. The dye labeled cDNA samples were hybridized to the whole human genome microarray (Cat. no. G4112F; Design ID: 014850, Agilent Technologies, USA). We followed the experimental instructions provided by the kit supplier for RNA amplification, cDNA preparation, dye labeling, and array hybridization.

Microarray Image Analyses and Data
Processing. The microarrays were scanned at 5 m resolution using a GenePix 4000B scanner (Molecular Devices, USA). The scanned images were analyzed with GenePix Pro 6.1, the software package provided by Molecular Devices, USA. Quality control of the images was performed with Acuity 4.0 (Molecular Devices, USA) software package. The following quality criteria have been used: (1) at least 90% of the pixels in the spot had intensity higher than background plus two standard deviations, (2) less than 2% saturated pixels in the spot, (3) signal to noise ratio (defined as ratio of the background subtracted mean pixel intensity to standard deviation of background) was 3 or above for each channel, (4) the spot diameter was between 45 and 70 m, (5) the regression coefficient of ratios of pixel intensity was 0.6 or above, and (6) feature has opposite change direction at dye-swap. Only features that passed all quality control criteria in all microarrays were analyzed further. To identify significantly expressed genes, the R package "limma" was used. "Within array" normalization was carried out with Lowess method, and "between arrays" normalization with "quantile" algorithm in the "limma" package. Features that did not pass the quality control criteria in all arrays were given weight of 0 and were excluded from the analyses. Differential expression and false discovery rate (FDR) were assessed using a linear model and empirical Bayes moderated F statistics. Genes that passed FDR ≤ 0.01 in all time points with a fold change of 2 or more in at least one time point and in the same direction in all time points were uploaded into MetaCore tool (Thomson Reuters, NY) for the identification of pathways affected by LE. The same data sets were uploaded to TM4, microarray data management, and analysis software to generate the heat map [4].

Results
The whole genome microarray (Agilent, USA) featured approximately 41,000, 60 nucleotide long oligonucleotide probes. These probes represent all known human transcripts. The microarray data were deposited in NCBI's "Gene Expression Omnibus" database repository; accession number is GSE44290 [11].
Seven genes differentially regulated in microarray were selected for validation using quantitative real time PCR. The direction of fold change in gene expression observed in microarray analyses matched with the real time PCR data ( Table 1).

Primary Screening of Microarray Data Yielded 17,135 Genes That Passed Technical QC in All Microarrays and Have Been
Analyzed Further. 17,135 features including replicated spots passed technical QC (see "Section 2.4") in all microarrays. Hybridization signals that passed all the quality control criteria as mentioned in the materials and methods were compared to the signals obtained by dye-swap hybridization. We combined biological replicates for dye swap in each treatment time point. By this method of data filtering, we selected 2,724, 2,730, and 5,463 genes for 24, 48, and 72 h, respectively (see Supplementary File S-1 in Supplementary Materials available online at http://dx.doi.org/10.1155/2013/471739). By this screening, we found that 1178 genes were common for all the three time points irrespective of the direction of expression change viz up or down regulation. Details of these common genes are listed in the Supplementary File S-1. In this data table, the same gene may be upregulated at one time point and downregulated at other time points. These 1,178 genes represented expression data that passed initial QC for all the three time points irrespective of direction of gene expression change triggered by LE treatment. Except for the QC described above, no preset fold change (FC) or false discovery rate (FDR) parameters were applied to these data sets. The lists of genes thus identified to be common for all three time points are detailed in "Spot test" sheet in the Supplementary File S-1. 24hFDR001, 48hFDR001, and 72hFDR001, respectively (Supplementary File S-1). The genes selected by this method were compared for commonly regulated genes between the different treatment time points. This analysis found 223, 191, and 2,009 genes to be unique to 24, 48, and 72 h, respectively. 38 genes were established to be common between 24 and 48 h, 489 between 24 and 72 hours, and 356 genes between 72 and 48 h. 221 genes were common for all the three time points. Though the genes may be common for the indicated parameters, the direction of regulation may not be the same for these genes at the compared time points. In this list, only the gene names were common. Other than FDR ≤ 0.01 the FC was not factored in this selection. Therefore, the same gene may have a positive and negative fold change at different time points (Sheet "FDR ≤ 0.01, " Supplementary File S-1). This data is graphically represented in Figure 1.

Application of Dual Filter to the Data Set Identified 82
Genes to Be Differentially Regulated by LE. Next, we limited the data analyses to genes that showed similar expression trend either up or down regulated over the three different experimental time points. This analysis yielded the following results. From across three different time points, 941 (516/425) out of 1,178 genes showed similar expression pattern (Sheet "FC Same trend, " Supplementary File S-1). From this list, 362 (248/114) genes showed FC ≥ ±2 (Sheet "FC ≥ ±2, " Supplementary File S-1). This list was further screened for genes with FDR ≤ 0.01 and thus we picked 89 (58/31) genes (Sheet "FDR ≤ 0.01 and FC ≥ ±2, " Supplementary File S-1). Of the 89 genes, 82 had known valid accession IDs. Thus, from 41,000 probe sets on Human Genome microarray, we successfully selected 82 (54/28) genes that showed similar trend across all the three time points with FDR 0.01. In this list, all the genes in at least one of the three points qualified FDR ≤ 0.01 and FC ≥ ±2 (Tables 2(a) and 2(b)). The relatedness of the expression of 89 genes that are significantly down (58) or up (31) regulated at all the three time points is depicted in the heat map ( Figure 2).

LE Regulates Key Pathways Involved in Immune
Response, Development, and So Forth. The 82 genes thus selected were uploaded to MetaCore tool in GeneGo suite (GeneGo, USA) to determine the possible signaling pathways affected by LE. This enrichment analysis identified 98 pathways to be regulated by LE. We grouped these 98 pathways that have similar effect. The top ten pathway maps of the 98 identified as above are given in Figure 3. A detailed list of all the 98 pathways affected by LE is given in Table 3. Overall, the data suggested LE predominantly affected immune response (#27), followed by metabolism (#17), developmental (#12), gene expression regulation (#6), cell adhesion (#5), cystic fibrosis transmembrane conductance regulation (CFTR), and chemotaxis (#4) pathways.
In this analysis, the transcription assembly of RNA Polymerase II preinitiation complex on TATA-less promoters pathway passed QC set at FDR ≤ 0.05. The value was 1.216 × 10 −5 ( Figure 4). Three objects from the uploaded data were found common among the 18 genes that are listed to be involved in this regulation. The key genes involved in this pathway are RNA pol II that is important for the transcription of TATA less promoters. This gene was underexpressed in LE treated cells. It suggests the deregulation of TATA less promoter genes. In genome wide analyses, only 15% of the human genes have been reported to be regulated by classical TATA box containing promoter. Redundancy/overlap of these regulations has been observed in several genes. Therefore, even though it appears that the regulation of TATA less gene transcription is altered, complete loss of their activity may be limited. Otherwise, LE treatment must result in unacceptably high toxicity, which was not encountered in animal studies [7]. Figure 4 illustrates the TF-II and the status of the three genes that are directly involved in this pathway.

Disease Biomarker Analysis Suggests LE Affects Vasculitis.
The expression of 82 genes that were identified to be influenced by LE was subjected to disease biomarker enrichment analysis in MetaCore tool. This investigation picked up a large number of diseases from this data set. The top 10 diseases along with the corresponding -values are listed in Table 4. It also contains the number of genes found in the microarray data that matched genes cataloged in Meta-Core disease biomarker enrichment database. It appears that LE predominantly affected genes/pathways involved in the inflammation of blood vessels, a condition called Vasculitis. Autoimmune disorders like arthritis and pancreatitis are attributed to be two of many etiology for vasculitis [12,13]. In this table, arthritis is the second most disease affected by LE ( ≤ 6.185 − 12). As reported in the literature, A. aspera is used as a cure for pancreatic cancer [7], as an antiarthritic medication [14], and to stop bleeding. The outcome of disease network analysis supports its medicinal use. It is notable that LE decreased the CD44 expression (Table 2(a)). The cell surface glycoprotein CD44 is involved in a wide variety of interactions that include receptors, cellcell interaction, adhesion, and migration (metastasis). Along with CD44, MMPs, VEGF, and hyaluronoglucosaminidase 1   Table 3.
(HYAL1; Table 2(a)), levels decreased in MiaPaCa-2 cells treated with LE. Of these molecules, the inhibition of HYAL1 enzyme was shown to decrease prostate cancer progression in vivo [15]. Therefore, LE is likely to suppresses tumor progression and metastasis as demonstrated earlier in our animal studies [7].

Discussion
In complementary and alternate medicine (CAM), unrefined (crude) preparations of natural products are used for the treatment of various illnesses. It is believed that compounds present in a crude preparation modulate the effect of the active molecule [16]. In this approach, the use of crude extract by CAM practitioners is akin to multidrug therapy. The difference is, in multidrug therapy the chemical nature as well the mechanisms of action of the component drugs are reasonably well studied. Even when combination chemotherapy is practiced, the added drugs are effective only in a few drug combinations [17]. Most of the modern day wonder drugs originated from or derivatives of the active molecules identified from natural products [18,19]. Once isolated, the active components from these natural products either end up being intolerably toxic, lose their activity, or needed in impractically large doses. Also, the drug discovery pathway is long, arduous, and expensive [20]. While it is ideal to identify the active molecule, we should note that a majority of the world population depends on natural remedies. Shortcomings in the use of natural products as medicine are the lack of quality control of these preparations, dosage optimization and contra indications. Various nodal agencies have embarked on an effort to validate the medicinal benefits, quality control, and dosage standardization of natural products that are widely used by CAM practitioners [1].
We are interested in the anticancer properties of A. aspera. It is noteworthy that LE modulates the expression of genes involved in embryonic development, cell adhesion, transcription factors that regulate cell growth, and so forth. These effects support its' use as abortifacient [21]. Hitherto, we have systematically demonstrated that A. aspera, a medicinal plant used in Ayurvedic medicine, inhibited the growth of cancer cells in vitro and orthotopic pancreatic tumor in vivo. Many parts of this plant are used as a cure for a wide array of ailments that include fertility control, fetal expulsion, wound healing, and cancer therapy. Since LE is a crude extract we expected a large number of genes would be affected in the treated cells. However, stringent data filtering determined just 89 genes were affected by LE. Of these, seven were probes used for quality control. Thus, effectively, we have narrowed down the number of genes affected by LE to 82. The pathway data analyses demonstrated predominantly immune response is affected by LE, followed by metabolism, development, transcription, and so forth. The effects on developmental, gene expression regulation, and cell adhesion pathways imply that LE directly affects cancer growth and metastasis.
Mutations in cystic fibrosis gene are identified as a risk factor for chronic pancreatitis and pancreatic ductal adenocarcinoma [22,23]. Molecular mechanistic study suggest loss of cystic fibrosis transmembrane conductance regulation (CFTR) and a corresponding increase in MUC4 expression is observed in about 81% of pancreatic ductal adenocarcinoma cell lines [24,25]. The microarray data suggest genes involved in CFTR regulation are modulated by LE. Though the effect of LE on CFTR pathway is highlighted in the gene expression studies, the molecular mechanism by which LE regulates CFTR or related gene function is not known. The short listing of CFTR and other top 10 pathways as affected by LE is worth  Figure 4: Canonical pathway analyses in MetaCore tool identified Transcription assembly of RNA polymerase II preinitiation complex on TATA-less promoters to be the top scored pathway map ( ≤ 1.204 − 05) in LE treated cells. Red thermometers show an object that is upregulated by LE. Blue thermometer show the objects downregulated by LE. The exact fold change obtained from three time points with the corresponding SD is given in Table 2. The big arrow indicates the "pathway start. " TR: transcriptional regulation; CS: complex subunit; B: binding; grey arrow: technical link; green arrows: positive effect blue arrows represent positive, red negative, and grey unspecified interactions. Boxes on lines denote the type of regulation where P is phosphorylation, B is binding, and TR is transcriptional regulation. A detailed legend for the objects in this figure is given in the Supplementary File SF-1. further investigations. The microarray data are consistent with our experimental findings that LE is preferentially cytotoxic to pancreatic cancer cell lines in vitro [8] and the use of LE in the treatment of pancreatic cancer by ayurvedic physicians.
Our results demonstrate that it is possible to pin point possible mechanism of action of crude extract using carefully planned experiments. We applied several assumptions in our data analyses. For instance, we eliminated expressions that did not follow the same pattern. There are instances where a group of genes may be up regulated in the short term and down regulated over long term exposure to LE or vice versa. Though less likely, another set of genes may be up and down regulated in a cyclical manner. By considering unidirectional regulation, we would have eliminated such unique changes effected by LE. Ambiguities in gene expression analysis may be reduced and the robustness of expression data may be increased by using several experimental replicates, multiple different cell lines, and different drug concentrations. Additionally, gene expression analyses may be performed in tumor samples obtained from animal model that had been treated with LE.
The identification of regulation of 82 genes from a large data set encourages us to think that natural products may target a limited set of genes. There are evidences to support such observations. For instance, almost identical results were obtained in MCF-7 breast cancer cells that were treated with a crude extract of Anoectochilus formosanus and a single compound drug plumbagin [26]. However, natural products may have much broader impact on global gene expression pattern. Therefore, we must be cautious in the universality of this approach.
In conclusion, we have demonstrated that LE affected the expression of only a few genes. Using similar methodology, it is possible to get reasonably good overview of the changes induced by crude preparation of natural products with medicinal properties.
Collaborative Translational and Clinical Studies component is supported by Grant no. 1UL1TR000460, the University of Miami CTSI, from the National Center for Advancing Translational Sciences and the National Institute on Minority Health and Health Disparities. Additional support was provided by Sheibar Foundation Grant to BA; NIH Grants: 1R01AT003544 and R01CA156776-01 to BLL.