Discovery of genetic biomarkers contributing to variation in drug response of cytidine analogues using human lymphoblastoid cell lines
BMC Genomics volume 15, Article number: 93 (2014)
Two cytidine analogues, gemcitabine and cytosine arabinoside (AraC), are widely used in the treatment of a variety of cancers with a large individual variation in response. To identify potential genetic biomarkers associated with response to these two drugs, we used a human lymphoblastoid cell line (LCL) model system with extensive genomic data, including 1.3 million SNPs and 54,000 basal expression probesets to perform genome-wide association studies (GWAS) with gemcitabine and AraC IC50 values.
We identified 11 and 27 SNP loci significantly associated with gemcitabine and AraC IC50 values, respectively. Eleven candidate genes were functionally validated using siRNA knockdown approach in multiple cancer cell lines. We also characterized the potential mechanisms of genes by determining their influence on the activity of 10 cancer-related signaling pathways using reporter gene assays. Most SNPs regulated gene expression in a trans manner, except 7 SNPs in the PIGB gene that were significantly associated with both the expression of PIGB and gemcitabine cytotoxicity.
These results suggest that genetic variation might contribute to drug response via either cis- or trans- regulation of gene expression. GWAS analysis followed by functional pharmacogenomics studies might help identify novel biomarkers contributing to variation in response to these two drugs and enhance our understanding of underlying mechanisms of drug action.
Both gemcitabine and AraC are widely used in the treatment of a variety of cancers and both display wide individual variation in drug response [1–6]. Pharmacogenomic studies have the potential to provide insight into mechanisms underlying individual variation in response to these two drugs [7–11]. Many previous pharmacogenetic studies focused on the bioactivation and metabolism pathways for cytidine analogues [12, 13]. For example, SNPs in genes encoding ribonucleotide reductase (RRM1) and cytidine deaminase (CDA) were found to be associated with gemcitabine chemosensitivity in the NCI-60 cell lines or with active gemcitabine metabolite plasma levels [14–16]. Those findings provided the initial evidence that genetic variation might contribute to variation in cytidine analogue response. We previously used the “Human Variation Panel”, a genomic data-rich lymphoblastoid cell line model system, to identify markers that might contribute to variation in response to these two cytidine analogues [17, 18]. These cell lines have proven to be a powerful tool for both the identification of pharmacogenomic hypotheses and for the pursuit of hypotheses from the clinical GWAS [19–21]. However, the earlier studies were performed with less dense SNP coverage, in the present study, we expanded our previous 550 K SNP data to include a total of 1.3 million SNPs obtained with both Illumina and Affymetrix SNP genotyping platforms in an attempt to identify additional genes or SNPs that might be associated with drug response. To follow-up the candidates, we performed functional studies using tumor cell lines in an attempt to determine possible underlying mechanisms that might help us to better understand mechanisms of action for these two drugs. The results of the comprehensive series of experiments described subsequently resulted in the identification of several novel SNPs and genes associated with gemcitabine and AraC drug response in these cell lines. These results could be tested in future clinical studies to determine whether they might help predict response to gemcitabine and AraC.
One hundred and seventy four human lymphoblastoid cell lines from 60 Caucasian-American (CA), 54 African-American (AA) and 60 Han Chinese-American (HCA) (sample sets HD100CAU, HD100AA, and HD100CHI) subjects were purchased from the Coriell Cell Repository (Camden, NJ). All of these cell lines had been obtained and anonymized by the National Institute of General Medical Sciences prior to deposit, and all subjects had provided written consent for the use of their DNA and cells for experimental purposes. Human SU86 pancreatic cancer cell lines were a gift from Dr. Daniel D. Billadeau (Department of Immunology and Division of Oncology Research, Mayo Clinic College of Medicine). Human breast cancer MDA-MB-231 and leukemia BDCM and THP-1 cell lines were purchased from the American Type Culture Collection (ATCC) (Manassas, VA) and were cultured in DMEM with 1% L-glutamine (Mediatech) supplemented with 10% FBS (Mediatech). Other cell lines were maintained in RPMI medium 1640 with 1% L-glutamine (Mediatech) supplemented with 10% FBS (Mediatech).
Drugs and cell proliferation assays
Gemcitabine was provided by Eli Lilly (Indianapolis, IN). AraC was purchased from Sigma-Aldrich (St. Louis, MO). Cytotoxicity assays with the lymphoblastoid and tumor cell lines were performed with the CellTiter 96® AQueous Non-Radioactive Cell Proliferation Assay (Promega Corporation, Madison, WI) as previously described . IC50 values were calculated using a three or four parameter logistic model with the R package “drc” (http://cran.r-project.org/web/packages/drc/drc.pdf), as described previously .
In order to validate the imputation results, six top imputed SNPs (rs10447475, rs4621668, rs11215400, rs10926784, rs3196512 and rs7762319) were selected for genotyping using Applied Biosystems TaqMan technology (Life Technologies, Grand Island, NY). One SNP (rs7762319) was not genotyped because the assay failed functional test, and four of the remaining five SNPs were successfully genotyped. Among these four SNPs, rs11215400 was a pre-designed assay, while the remaining three SNPs (rs10447475, rs4621668 and rs10926784) were customized assays designed with Custom TaqMan® Assay Design Tool (Life Technologies, Grand Island, NY). Primer and probe sequences for these assays are available upon request. PCR protocols were followed according to the manufacturer’s guidelines for the 384-well format. PCR amplifications were performed using Applied Biosystems® TaqMan® Genotyping Master Mix with an Applied Biosystems® Veriti® 384-Well Thermal Cycler (Life Technologies, Grand Island, NY), and PCR products were analyzed on an Applied Biosystems 7900HT .
Transient transfection and RNA interference
Specific siGENOME siRNA SMARTpool® reagents against a given gene, as well as a negative control, siGENOME Non-Targeting siRNA Pool #2, were purchased from Dharmacon Inc. (Lafayette, CO). Human pancreatic cancer SU86 and breast cancer MDA-MB-231 cell lines were used to perform the siRNA knockdown studies. The lipofectamine RNAiMAX transfection reagent (Invitrogen, Carlsbad, CA) was used for siRNA reverse or forward transfection. Specifically, cells were seeded into 96-well plates and were mixed with siRNA-complex consisting of 20–50 nM of specific siGENOME siRNA SMARTpool or non-targeting negative control (Dharmacon) and the lipofectamine RNAiMAX transfection reagent. The human leukemia cell lines, BDCM and THP-1, were transfected with electroporation using the Nucleofector System with 500 nM of specific or negative siRNA (Lonza Inc., Basel, Switzerland).
Quantitative reverse transcription-PCR (QRT-PCR)
Total RNA was isolated from cultured cells with the Qiagen RNeasy kit (QIAGEN Inc. Valencia, CA), followed by QRT-PCR performed with the 1-step, Brilliant SYBR Green QRT-PCR master mix kit (Stratagene, La Jolla, CA). Specifically, primers purchased from Qiagen were used to perform QRT-PCR with the Stratagene Mx3005P™ Real-Time PCR detection system (Stratagene). All experiments were performed in triplicate with β-actin as an internal control. Reverse transcribed Universal Human reference RNA (Stratagene) was used to generate a standard curve. Control reactions lacked RNA template.
Caspase-3/7 activity assay
Caspase-3/7 activity was measured with the Caspase-Glo®3/7 Assay kit (Promega). Specifically, siRNA-transfected cells (100 μl) were seeded overnight into 96-well plates at a density of 10,000 cells per well and were then treated with DMSO or increasing concentrations of gemcitabine or AraC for 48 h. 100 μL of Caspase-Glo® 3/7 Reagent was then added to each well, and the cells were incubated at room temperature for 1 h, followed by the measurement of luminescence with a Safire2 microplate reader (Tecan Trading AG, Switzerland). The luminescent signal was proportional to caspase-3/7 activity and was used as a measure of apoptosis. Wells containing only culture medium were used as controls.
Cancer cignal finder 10-pathway reporter array
The Cignal Finder Arrays consist of 10 dual-luciferase reporter assays for ten cancer-related signaling pathways. Each reporter construct is a mixture of an inducible transcription factor (TF) responsive firefly luciferase reporter and a constitutively expressing Renilla construct at a ratio of 40:1, respectively (SABioscience Co., Frederick, MD). Specifically, cells were reversely transfected with 30 nM of specific siRNA pools in 96-well plates using Lipofectamine RNAiMAX reagent (Invitrogen) for 24 h, followed by transfection with 100 ng of each reporter construct. Forty-eight h after the transfection, a dual-luciferase assay was performed with the Dual-Luciferase Reporter Assay System (Promega) in a Safire2 microplate reader (Tecan).
Electrophoresis mobility shift assays (EMSA)
Based on the genome-wide association results, we performed EMSA for SNPs in potential regulatory regions of genes that were associated with the measured phenotypes. Specifically, double-strand probes were 5′-end labeled with biotin and electrophoresis was performed with 5% acrylamide gels, followed by autoradiography. Competition experiments were performed with excess non-labeled probe.
Genome-wide gene expression and SNP analysis
Expression array data were obtained for all 174 lymphoblastoid cell lines (LCLs) as previously described . Illumina HumanHap550K and 510S BeadChips, which assayed 561,298 and 493,750 SNPs, respectively, were used to obtain genome-wide SNP data for these LCLs . Genotyping was performed in the Genotype Shared Resource (GSR) at the Mayo Clinic, Rochester, MN. We also obtained publicly available Affymetrix SNP Array 6.0 Chip SNP data which involved 643,600 SNPs unique to the Affymetrix SNP array for the same cell lines. After quality control (QC), SNPs with call rates <0.95, Hardy-Weinberg Equilibrium (HWE) P values < 0.001, or MAFs <5% were excluded, as were DNA samples with call rates <0.95. A total of 1,348,798 SNPs that passed QC were used to perform the association studies.
SNPs not genotyped were imputed across a region 200 kb up or downstream of the selected genes harboring or close to the SNPs associated with drug response in the LCLs. Imputation was performed using Beagle (v3.3.1)  with the 11/23/2010 release of the 1000 Genomes project as a reference population . Imputed SNPs with a dosage R2 quality measure of less than 0.3, and SNPs with MAF <0.01 were not included in the analysis. Four of the imputed SNPs were genotyped for validation, the average squared difference between the count of the same allele in the imputed and genotyped version of these SNPs was computed to measure the concordance of the imputed genotype with actual genotype, a smaller difference indicating greater concordance.
Partial Pearson correlations were used to quantify the association between: SNPs and mRNA expression; SNPs and IC50; and mRNA expression and IC50. IC50 was transformed to remove skewness using a log transformation for gemcitabine and van der Waerden rank transformation for AraC. The adjustment variables in the partial correlation were race and gender if SNPs were not involved; or race, gender and five eigenvectors controlling for population stratification as described previously . These partial correlations were tested using a Wald test, false discovery q-values  were also computed for each test.
Genome-wide SNP vs. drug cytotoxicity association study and imputation analysis
Previously, we had performed GWAS using only the 550 k SNP data set for this cell line system . In the current study, we expand the SNPs studies to include additional Illumina SNPs as well as publically available SNP data obtained with Affymetrix 6.0 SNP data for the same cell lines to identify additional novel potential pharmacogenomic biomarkers. As a result, we performed an analysis for the association of 1,348,798 SNPs with IC50 values for gemcitabine and AraC (Figure 1A and 1B). The most significant SNP for gemcitabine was rs1598848 with a P value 7.08 × 10-7 (r = −0.391, MAF = 0.473), while the most significant SNP for AraC was rs4078252 with a P value 1.54 × 10-7 (r = 0.405, MAF = 0.198) (Additional file 1: Table S1A and B). Fourteen SNPs for gemcitabine and 34 for AraC had P values <10-5, and 143 SNPs for gemcitabine and 204 SNPs for AraC had P values <10-4, respectively. One hundred and twenty six SNPs with P < 10-3 were common to both drugs (Additional file 1: Table S1C). To explore ungenotyped SNPs that might be functional, we imputed SNPs surrounding the selected genes (+/−200 kb) harboring or close to the most significant SNP loci using 1000 Genomes data as a reference (Additional file 1: Figure S1A, B, and C). Besides the “observed SNPs” on the genotyping platforms, there were 23 imputed SNPs for gemcitabine and 35 for AraC, respectively, that were also associated with drug response IC50 values (P < 10-4) (Additional file 1: Table S2A and B). In order to determine the accuracy of imputation, we selected 6 imputed SNPs (rs10447475, rs4621668, rs11215400, rs10926784, rs3196512 and rs7762319) that were among the top two SNPs from each gene region associated with drug response IC50s (P < 10-3) to genotype using Taqman assay. Four SNPs ((rs10447475, rs4621668, rs11215400 and rs10926784) were successfully genotyped. The average squared difference between the count of the same allele in the imputed and genotyped version of these 4 SNPs was low ranging from 0.02-0.065 indicating that the concordance was high (Additional file 1: Figure S2).
“Integrated analyses” of SNP loci, basal expression and drug cytotoxicity
We also performed integrated analyses of SNPs, expression array and cytotoxicity data [18, 23]. To do that, we began with SNPs that had P values <10-3. We selected a less stringent P value cutoff to include as many potential candidates as possible for follow-up functional genomic studies. Next, we tested expression probe sets that were associated with these SNPs, followed by association of those expression probe sets with drug IC50 values, ie we performed an “integrated analysis”. In these analyses, we used SNP loci, defined as a region that contained at least 2 SNPs with P < 10-4 or 1 SNP with P < 10-4, plus 3 additional SNPs with P < 10-3 within +/−100 kb surrounding the most significant SNPs. All of the SNPs within each of those loci are listed in Additional file 1: Table S2, which includes genotyped as well as imputed SNPs. We identified 11 loci containing 166 SNPs for gemcitabine and 23 loci with 187 SNPs for AraC that were associated with IC50 values for these two drugs, respectively (Table 1). Four loci containing 4 genes – HLA-DRA, ZNF215, MASS1, and PLD5 – were common to both drugs (Table 1).
The integrated analyses identified 66 SNPs in 6 loci that were associated with gemcitabine IC50 values and the expression of 12 genes represented by 17 probesets. Those 17 probesets were also associated with gemcitabine IC50 values (P < 10-4). The integrated analyses also identified 36 SNPs in 3 loci that were associated with AraC IC50 values and the expression of 9 genes (10 probesets) (Table 2). For gemcitabine, 19 SNPs were within cis-regulatory regions for PIGB or C3orf23. No cis- regulation between SNP and gene expression was identified for AraC. Of interest, SNPs in PIGB were associated with the expression of that gene (lowest P = 5.97 × 10-9) as well as the expression of FKBP5 (lowest P = 2.70 × 10-6), a gene that we previously reported to play an important role in response to gemcitabine and AraC as well as many other chemotherapeutic agents including gemcitabine and AraC [17, 27]. We next moved to further analyses of candidate genes identified during the integrated analysis.
Follow-up functional validation of candidate genes in cancer cells
Since the regulation of gene expression is tissue specific , we wanted to functionally validate in cancer cell lines candidate genes selected based on our analysis performed with LCLs. The tumor cell lines that we selected were based on the expression of the genes of interest and on the clinical uses of these two drugs. Gemcitabine is used in the treatment of pancreatic cancer but it is also used to treat other solid tumors such as breast cancer, while AraC is first-line therapy for acute myelogenous leukemia (AML). Therefore, we selected one human pancreatic cancer cell line, SU86, one breast cancer cell line, MDA-MB-231 and two leukemia cell lines, BDCM and THP1, to functionally characterize the genes of interest. Twenty-six genes were selected based on a series of criteria including association P value, SNP locus, whether the gene was expressed in LCLs and the biological function of the genes (Table 3 and Figure 2). To determine the functional impact of those genes, we used specific siRNA pools to knockdown the 26 candidate genes, followed by QRT-PCR and MTS cytotoxicity assay. Eleven genes showed an effect on gemcitabine cytotoxicity, 10 on AraC and 5 were common to both drugs. Knockdown of PIGB, ZADH2, PSME3, DOK6, TGFBI, and HLA-DRA in both SU86 and MDA-MB-231 cells significantly desensitized the cells to gemcitabine (Table 3 and Figure 3), consistent with the association study results. Knockdown of C4orf169, TUSC3, LNX2, RIS1, and SMC2 and HLA-DRA in both SU86 and MDA-MB-231 cells significantly desensitized the cells to AraC (Table 3 and Figure 4). Finally, knockdown of HLA-DRA in THP-1 leukemia cells, LNX2 in BDCM cells, and SMC2 and RIS1 in both THP1 and BDCM cells also desensitized the cells to AraC, results that were also consistent with our association results (Figure 4).
We next wanted to determine whether the cytotoxic effects of those genes might involve apoptosis. Therefore, we performed caspase-3/7 activity assays after knockdown of the candidate genes in SU86 cells. As shown in Figure 5A and 5B, down-regulation of PIGB, DOK6, TGFBI, ZADH2, PSME3, and HLA-DRA in SU86 cells significantly decreased caspase-3/7 activity after treatment with gemcitabine as compared with negative control siRNA-treated cells. Similar results were also observed for AraC treatment following siRNA knockdown of TUSC3, C14orf169, and HLA-DRA. However, knockdown of LNX2, RIS1, and SMC2 did not alter the cellular caspase-3/7 activity (Table 4), suggesting that a different mechanism was involved.
Finally, we used the Cancer Cignal Finder Array (SABioscience) that consists of 10 dual-luciferase reporter gene assays to determine whether our candidate genes might affect any of the 10 cancer-related signaling pathways in SU86 cells by measuring changes in transcriptional activities of 10 key transcription factors (TF) after knockdown of each candidate gene. We observed changes in transcriptional activity of several TFs after knockdown of specific genes in SU86 cells, suggesting that these genes might be involved in the regulation of a particular cancer-related signaling pathway or pathways that might contribute to resistance to gemcitabine and AraC (Figure 6 and Table 3). For example, knockdown of PIGB resulted in a decrease in transcriptional activity of Elk-1/SRF, AP1, NFκB, and Myc/MAX in SU86 cells, indicating a down-regulation of these signaling pathways. Knockdown of DOK6 dramatically decreased the transcription activities of both NFκB and AP1 in the NFκB and MAPK/JNK pathways, while the activity of the transcription factor Myc/MAX that is involved in the c-Myc pathway was increased significantly after ZADH2 knockdown. However, we did not observe any significant changes after SMC2 knockdown.
Functional characterization of PIGB SNPs
When we performed integrated analysis among SNPs, gene expression and gemcitabine cytotoxicity, we found that the only cis- regulated SNPs mapped to PIGB. Knockdown of PIGB resulted in desensitization of cancer cells to gemcitabine. PIGB contained 7 SNPs that were associated both with gemcitabine response (P < 10-3) and with its own gene expression (P < 10-4) (Table 5). PIGB expression was also significantly correlated with gemcitabine cytotoxicity (P = 8.95 × 10-5 and P = 5.31 × 10-3 for two different probe sets for PIGB mRNA). We also determined LD patterns for those 7 SNPs using HapMap data for each ethnic group. As shown in Figure 7A, LD patterns differed among the three ethnic groups. In both CHB/JPT and CEPH groups, those 7 SNPs were in tight LD, while there was not significant linkage among the SNPs in the YRI population. The top 3 SNPs in PIGB, including rs2290344, a nonsynonymous coding SNP (M161T) in exon 4, rs28668016 in the 5′-UTR, and rs11636687 in the 5′ flanking region (Table 5) were selected for further functional characterization.
We first determined PIGB expression levels in 37 LCLs selected on the basis of genotypes for those 3 SNPs using both QRT-PCR assay and expression array data to confirm the association between the SNPs and PIGB expression. Cells carrying the variant alleles showed significantly lower expression levels than did WT cells (Figure 7B). We next determined the functional impact of these 3 SNPs. As shown in Figure 7C, overexpression of a construct for the PIGB coding SNP (rs2290344, M161T) in SU86 and MDA-MB-231 cells did not alter either mRNA or protein levels, nor did it have an effect on gemcitabine cytotoxicity. We then determined whether the two SNPs in regulatory regions, rs11636687 and rs28668016, might have functional impact. We performed electrophoresis mobility shift assays (EMSAs) for these two SNPs to determine whether there might be differences in binding patterns for possible transcription factors. Interestingly, the results from EMSA showed that DNA-protein binding was significantly increased for the probe containing the variant sequences for these two SNPs in both SU86 and MDA-MB-231 cells (Figure 7D). These results suggested that these two SNPs might alter the binding of transcription factors and, as a result, affect PIGB expression level.
We previously performed a genome wide SNP association study with 550 K SNPs obtained with Illumina HumanHap550 BeadChips for the same cell lines to identify common polymorphisms that might influence both gene expression and response to these two drugs . In the present study, we expanded the number of SNPs from 550 K to include over 1.3 million SNPs and selected candidate genes for functional follow-up studies based on SNP loci. This dense SNP coverage made it possible to identify many more candidates for functional follow-ups. That enabled us to take a different approach by focusing on “SNP loci” instead of single SNPs. The results listed in Table 3 show that 11 of 26 candidate genes selected in this fashion were validated functionally, while only two other genes from the previous 550 k studies were functionally validated .
We also tested the concordance of the results generated with 550 K and 1.3 million SNPs if we had used the same strategy as we did in the current study, i.e. using SNP loci to perform the association studies. The majority of top SNP peaks from the 550 K SNP data for both drugs displayed less significant SNPs for each locus as compared to the 1.3 million SNP data (Additional file 1: Table S3). These observations illustrate the advantage of the present selection strategy for candidate identification, as well as the advantage of using denser SNP coverage.
Of the 26 candidate genes that we identified for further siRNA screening followed by MTS assay, eleven candidate genes, including PIGB, TGFBI, DOK6, PSME3, ZADH2, TUSC3, C14orf169, SMC2, LNX2, RIS1, and HLA-DRA, showed a significant effect on response to gemcitabine and/or AraC in SU86 and/or MDA-MB231 cells. To identify potential pathways with which these genes might be involved, we used a dual luciferase reporter gene assay to assess the impact of these genes on 10 major cancer-related signaling pathways. As shown in Figure 6 and Table 3, except for the SMC2 gene, knockdown of the other 10 genes in SU86 cells significantly altered activities, based on the luciferase assay for at least one of the 10 cancer related signaling pathways. Genes such as TGFB1 showed changes for the most pathways. While TGFB1 has been well studied, genes such as C14orf169, an unknown gene, also showed increased activity in 7 of the 10 pathways.
We also observed that the activities of the Elk-1/SRF, AP1, NFκB, and Myc/MAX pathways were significantly decreased in SU86 cells when PIGB was down-regulated by a specific siRNA. PIGB, a gene of the phosphatidylinositol glycan (PIG) class B, encodes an enzyme involved in the synthesis of a glycosylphosphatidylinositol (GPI) anchor that is a membrane attachment structure for many proteins, including membranous enzymes, receptors, differentiation antigens, and other biologically active proteins . GPI anchoring is essential for the expression of many of those proteins in either biological processes or cancer progression [30, 31]. The PIGB protein is a GPI mannosyltransferase III and is required for the transfer of the third mannose into the core structure of the GPI anchor [29, 32]. Previous studies have demonstrated that other PIG class members, such as PIGU and PIGT, are oncogenes in either human bladder cancer or breast cancer, respectively [33, 34]. Our findings indicate that PIGB is involved in sensitizing cancer cells to both gemcitabine and AraC, suggesting a possible role in oncogenic pathways as well as chemoresistance. The 8 PIGB SNPs were also associated with the expression of FKBP5, a gene that we previously reported to be important for gemcitabine and AraC response [17, 27]. Furthermore, PIGB expression itself is also correlated with FKBP5 gene expression. Although down regulation of PIGB altered FKBP5 mRNA level, overexpression of FKBP5 in PIGB stable knockdown cell lines did not change response to gemcitabine or AraC (Additional file 1: Figure S3). These observations indicate that PIGB influences the cytotoxicity of the two cytidine analogues through mechanisms that differ from FKBP5, in spite of the correlation of their expression levels observed in the LCLs. The exact mechanisms by which PIGB affects gemcitabine and AraC cytotoxicity need to be explored in the course of future experiments.
In addition to the characterization of candidate genes, we also focused on SNPs in the PIGB gene that showed cis-regulation of PIGB expression. SNPs in regulatory regions can influence drug response through an influence on gene expression. During our analysis, we found that most SNP associations with expression were through trans- regulation. The reason that we focused on SNPs in PIGB is because those SNPs displayed cis- regulations of PIGB and knockdown of PIGB showed an effect on cytotoxicity. The EMSA results also demonstrated “shifts” for the variant SNP sequences (Figure 7D), suggesting that PIGB gene expression might be regulated through binding to those transcription factors.
Previous studies demonstrated that one mechanism by which SNPs might influence drug cytotoxicity is through transcription regulation in either a cis- or trans-manner [18, 35–37]. In this analysis, we found SNPs that could both have cis or trans relationship. In addition to the SNPs that cis regulate PIGB, we also found that SNPs close to C3orf23 were not only cis-associated with its own gene expression, but also trans-correlated with the expression of ZADH2 which was confirmed to affect drug response of gemcitabine in our functional validation study. How those genetic variations located in the upstream of C3orf23 affect the expression of ZADH2 gene in a trans- manner remains unknown. One mechanism might be that those SNPs nearby C3orf23 could alter DNA sequence binding to transcription factors (TFs), microRNA, or other long non-coding RNA (lnc RNA), thus affect transcriptional regulation of their target genes including ZADH2 gene, which could in turn, affect gemcitabine response.
In summary, this study performed with LCLs followed by functional characterization has enhanced our understanding of the action of gemcitabine and AraC in the therapy of cancer. Although there are limitations associated with the use of LCLs [38, 39], this system has proven to be extremely useful, both to generate pharmacogenomic hypothesis and to test pharmacogenomic signals identified during the clinical GWAS [19–21]. Future studies using patient samples will now be required to confirm the candidates identified during this study.
Lymphoblastoid cell line
Electrophoresis mobility shift assays
Quantitative reverse transcription-PCR
Acute myelogenous leukemia
Genome-wide association studies
Tomasson MH, Xiang Z, Walgren R, Zhao Y, Kasai Y, Miner T, Ries RE, Lubman O, Fremont DH, McLellan MD, Payton JE, Westervelt P, DiPersio JF, Link DC, Walter MJ, Graubert TA, Watson M, Baty J, Heath S, Shannon WD, Nagarajan R, Bloomfield CD, Mardis ER, Wilson RK, Ley TJ: Somatic mutations and germline sequence variants in the expressed tyrosine kinase genes of patients with de novo acute myeloid leukemia. Blood. 2008, 111: 4797-4808. 10.1182/blood-2007-09-113027.
Heinemann V, Xu YZ, Chubb S, Sen A, Hertel LW, Grindey GB, Plunkett W: Cellular elimination of 2′,2′-difluorodeoxycytidine 5′-triphosphate: a mechanism of self-potentiation. Cancer Res. 1992, 52: 533-539.
Heinemann V, Hertel LW, Grindey GB, Plunkett W: Comparison of the cellular pharmacokinetics and toxicity of 2′,2′-difluorodeoxycytidine and 1-beta-D-arabinofuranosylcytosine. Cancer Res. 1988, 48: 4024-4031.
Kindler HL: In focus: advanced pancreatic cancer. Clin Adv Hematol Oncol. 2005, 3: 420-422.
Kern W, Estey EH: High-dose cytosine arabinoside in the treatment of acute myeloid leukemia: Review of three randomized trials. Cancer. 2006, 107: 116-124. 10.1002/cncr.21543.
Hertel LW, Boder GB, Kroin JS, Rinzel SM, Poore GA, Todd GC, Grindey GB: Evaluation of the antitumor activity of gemcitabine (2′,2′-difluoro-2′-deoxycytidine). Cancer Res. 1990, 50: 4417-4422.
Sabeti PC, Varilly P, Fry B, Lohmueller J, Hostetter E, Cotsapas C, Xie X, Byrne EH, McCarroll SA, Gaudet R, Schaffner SF, Lander ES, Frazer KA, Ballinger DG, Cox DR, Hinds DA, Stuve LL, Gibbs RA, Belmont JW, Boudreau A, Hardenbol P, Leal SM, Pasternak S, Wheeler DA, Willis TD, Yu F, Yang H, Zeng C, Gao Y, Hu H, et al: Genome-wide detection and characterization of positive selection in human populations. Nature. 2007, 449: 913-918. 10.1038/nature06250.
Altshuler D, Daly MJ, Lander ES: Genetic mapping in human disease. Science. 2008, 322: 881-888. 10.1126/science.1156409.
Baum AE, Akula N, Cabanero M, Cardona I, Corona W, Klemens B, Schulze TG, Cichon S, Rietschel M, Nothen MM, Georgi A, Schumacher J, Schwarz M, Abou Jamra R, Hofels S, Propping P, Satagopan J, Detera-Wadleigh SD, Hardy J, McMahon FJ: A genome-wide association study implicates diacylglycerol kinase eta (DGKH) and several other genes in the etiology of bipolar disorder. Mol Psychiatry. 2008, 13: 197-207. 10.1038/sj.mp.4002012.
Wang L, McLeod HL, Weinshilboum RM: Genomics and drug response. N Engl J Med. 2011, 364: 1144-1153. 10.1056/NEJMra1010600.
Wang L, Weinshilboum RM: Pharmacogenomics: candidate gene identification, functional validation and mechanisms. Hum Mol Genet. 2008, 17: R174-R179. 10.1093/hmg/ddn270.
Kim SR, Saito Y, Maekawa K, Sugiyama E, Kaniwa N, Ueno H, Okusaka T, Morizane C, Yamamoto N, Ikeda M, Yoshida T, Minami H, Furuse J, Ishii H, Saijo N, Kamatani N, Ozawa S, Sawada J: Thirty novel genetic variations in the SLC29A1 gene encoding human equilibrative nucleoside transporter 1 (hENT1). Drug Metab Pharmacokinet. 2006, 21: 248-256. 10.2133/dmpk.21.248.
Kim SR, Saito Y, Maekawa K, Sugiyama E, Kaniwa N, Ueno H, Okusaka T, Ikeda M, Morizane C, Yamamoto N, Yoshida T, Kamatani N, Furuse J, Ishii H, Saijo N, Ozawa S, Sawada J: Twenty novel genetic variations and haplotype structures of the DCK gene encoding human deoxycytidine kinase (dCK). Drug Metab Pharmacokinet. 2008, 23: 379-384. 10.2133/dmpk.23.379.
Kwon WS, Rha SY, Choi YH, Lee JO, Park KH, Jung JJ, Kim TS, Jeung HC, Chung HC: Ribonucleotide reductase M1 (RRM1) 2464G>A polymorphism shows an association with gemcitabine chemosensitivity in cancer cell lines. Pharmacogenet Genomics. 2006, 16: 429-438. 10.1097/01.fpc.0000204999.29924.da.
Sugiyama E, Kaniwa N, Kim SR, Kikura-Hanajiri R, Hasegawa R, Maekawa K, Saito Y, Ozawa S, Sawada J, Kamatani N, Furuse J, Ishii H, Yoshida T, Ueno H, Okusaka T, Saijo N: Pharmacokinetics of gemcitabine in Japanese cancer patients: the impact of a cytidine deaminase polymorphism. J Clin Oncol. 2007, 25: 32-42.
Yonemori K, Ueno H, Okusaka T, Yamamoto N, Ikeda M, Saijo N, Yoshida T, Ishii H, Furuse J, Sugiyama E, Kim SR, Kikura-Hanajiri R, Hasegawa R, Saito Y, Ozawa S, Kaniwa N, Sawada J: Severe drug toxicity associated with a single-nucleotide polymorphism of the cytidine deaminase gene in a Japanese cancer patient treated with gemcitabine plus cisplatin. Clin Cancer Res. 2005, 11: 2620-2624. 10.1158/1078-0432.CCR-04-1497.
Li L, Fridley B, Kalari K, Jenkins G, Batzler A, Safgren S, Hildebrandt M, Ames M, Schaid D, Wang L: Gemcitabine and cytosine arabinoside cytotoxicity: association with lymphoblastoid cell expression. Cancer Res. 2008, 68: 7050-7058. 10.1158/0008-5472.CAN-08-0405.
Li L, Fridley BL, Kalari K, Jenkins G, Batzler A, Weinshilboum RM, Wang L: Gemcitabine and arabinosylcytosin pharmacogenomics: genome-wide association and drug response biomarkers. PLoS One. 2009, 4: e7765-10.1371/journal.pone.0007765.
Ingle JN, Schaid DJ, Goss PE, Liu M, Mushiroda T, Chapman JA, Kubo M, Jenkins GD, Batzler A, Shepherd L, Pater J, Wang L, Ellis MJ, Stearns V, Rohrer DC, Goetz MP, Pritchard KI, Flockhart DA, Nakamura Y, Weinshilboum RM: Genome-wide associations and functional genomic studies of musculoskeletal adverse events in women receiving aromatase inhibitors. J Clin Oncol. 2010, 28: 4674-4682. 10.1200/JCO.2010.28.5064.
Liu M, Ingle JN, Fridley BL, Buzdar AU, Robson ME, Kubo M, Wang L, Batzler A, Jenkins GD, Pietrzak TL, Carlson EE, Goetz MP, Northfelt DW, Perez EA, Williard CV, Schaid DJ, Nakamura Y, Weinshilboum RM: TSPYL5 SNPs: association with plasma estradiol concentrations and aromatase expression. Mol Endocrinol. 2013, 27: 657-670. 10.1210/me.2012-1397.
Liu M, Wang L, Bongartz T, Hawse JR, Markovic SN, Schaid DJ, Mushiroda T, Kubo M, Nakamura Y, Kamatani N, Goss PE, Ingle JN, Weinshilboum RM: Aromatase inhibitors, estrogens and musculoskeletal pain: estrogen-dependent T-cell leukemia 1A (TCL1A) gene-mediated regulation of cytokine expression. Breast Cancer Res. 2012, 14: R41-10.1186/bcr3137.
Abo R, Hebbring S, Ji Y, Zhu H, Zeng ZB, Batzler A, Jenkins GD, Biernacka J, Snyder K, Drews M, Fiehn O, Fridley B, Schaid D, Kamatani N, Nakamura Y, Kubo M, Mushiroda T, Kaddurah-Daouk R, Mrazek DA, Weinshilboum RM: Merging pharmacometabolomics with pharmacogenomics using ‘1000 Genomes’ single-nucleotide polymorphism imputation: selective serotonin reuptake inhibitor response pharmacogenomics. Pharmacogenet Genomics. 2012, 22: 247-253. 10.1097/FPC.0b013e32835001c9.
Niu N, Qin Y, Fridley BL, Hou J, Kalari KR, Zhu M, Wu TY, Jenkins GD, Batzler A, Wang L: Radiation pharmacogenomics: a genome-wide association approach to identify radiation response biomarkers using human lymphoblastoid cell lines. Genome Res. 2010, 20: 1482-1492. 10.1101/gr.107672.110.
Browning SR, Browning BL: Rapid and accurate haplotype phasing and missing-data inference for whole-genome association studies by use of localized haplotype clustering. Am J Hum Genet. 2007, 81: 1084-1097. 10.1086/521987.
Abecasis GR, Auton A, Brooks LD, DePristo MA, Durbin RM, Handsaker RE, Kang HM, Marth GT, McVean GA, 1000 Genomes Project Consortium: An integrated map of genetic variation from 1,092 human genomes. Nature. 2012, 491: 56-65. 10.1038/nature11632.
Storey JD, Tibshirani R: Statistical significance for genomewide studies. Proc Natl Acad Sci U S A. 2003, 100: 9440-9445. 10.1073/pnas.1530509100.
Pei H, Li L, Fridley BL, Jenkins GD, Kalari KR, Lingle W, Petersen G, Lou Z, Wang L: FKBP51 affects cancer cell response to chemotherapy by negatively regulating Akt. Cancer Cell. 2009, 16: 259-266. 10.1016/j.ccr.2009.07.016.
Dimas AS, Deutsch S, Stranger BE, Montgomery SB, Borel C, Attar-Cohen H, Ingle C, Beazley C, Gutierrez Arcelus M, Sekowska M, Gagnebin M, Nisbett J, Deloukas P, Dermitzakis ET, Antonarakis SE: Common regulatory variation impacts gene expression in a cell type-dependent manner. Science. 2009, 325: 1246-1250. 10.1126/science.1174148.
Takahashi M, Inoue N, Ohishi K, Maeda Y, Nakamura N, Endo Y, Fujita T, Takeda J, Kinoshita T: PIG-B, a membrane protein of the endoplasmic reticulum with a large lumenal domain, is involved in transferring the third mannose of the GPI anchor. EMBO J. 1996, 15: 4254-4261.
Masterson WJ, Doering TL, Hart GW, Englund PT: A novel pathway for glycan assembly: biosynthesis of the glycosyl-phosphatidylinositol anchor of the trypanosome variant surface glycoprotein. Cell. 1989, 56: 793-800. 10.1016/0092-8674(89)90684-3.
Doering TL, Masterson WJ, Hart GW, Englund PT: Biosynthesis of glycosyl phosphatidylinositol membrane anchors. J Biol Chem. 1990, 265: 611-614.
Sutterlin C, Escribano MV, Gerold P, Maeda Y, Mazon MJ, Kinoshita T, Schwarz RT, Riezman H: Saccharomyces cerevisiae GPI10, the functional homologue of human PIG-B, is required for glycosylphosphatidylinositol-anchor synthesis. Biochem J. 1998, 332 (Pt 1): 153-159.
Guo Z, Linn JF, Wu G, Anzick SL, Eisenberger CF, Halachmi S, Cohen Y, Fomenkov A, Hoque MO, Okami K, Steiner G, Engles JM, Osada M, Moon C, Ratovitski E, Trent JM, Meltzer PS, Westra WH, Kiemeney LA, Schoenberg MP, Sidransky D, Trink B: CDC91L1 (PIG-U) is a newly discovered oncogene in human bladder cancer. Nat Med. 2004, 10: 374-381. 10.1038/nm1010.
Wu G, Guo Z, Chatterjee A, Huang X, Rubin E, Wu F, Mambo E, Chang X, Osada M, Sook Kim M, Moon C, Califano JA, Ratovitski EA, Gollin SM, Sukumar S, Sidransky D, Trink B: Overexpression of glycosylphosphatidylinositol (GPI) transamidase subunits phosphatidylinositol glycan class T and/or GPI anchor attachment 1 induces tumorigenesis and contributes to invasion in human breast cancer. Cancer Res. 2006, 66: 9829-9836. 10.1158/0008-5472.CAN-06-0506.
Kim DH, Jung HD, Kim JG, Lee J-J, Yang D-H, Park YH, Do YR, Shin HJ, Kim MK, Hyun MS, Sohn SK: FCGR3A gene polymorphisms may correlate with response to frontline R-CHOP therapy for diffuse large B-cell lymphoma. Blood. 2006, 108: 2720-2725. 10.1182/blood-2006-01-009480.
Ryan BM, Robles AI, Harris CC: Genetic variation in microRNA networks: the implications for cancer research. Nat Rev Cancer. 2010, 10: 389-402. 10.1038/nrc2867.
Wapinski O, Chang HY: Long noncoding RNAs and human disease. Trends Cell Biol. 2011, 21: 354-361. 10.1016/j.tcb.2011.04.001.
Feng WH, Hong G, Delecluse HJ, Kenney SC: Lytic induction therapy for Epstein-Barr virus-positive B-cell lymphomas. J Virol. 2004, 78: 1893-1902. 10.1128/JVI.78.4.1893-1902.2004.
Feng WH, Israel B, Raab-Traub N, Busson P, Kenney SC: Chemotherapy induces lytic EBV replication and confers ganciclovir susceptibility to EBV-positive epithelial cell tumors. Cancer Res. 2002, 62: 1920-1926.
This study was supported in part by U.S. National Institutes of Health research grants K22 CA130828, R01 CA138461, P50 CA102701, U19 GM61388 (The Pharmacogenomics Research Network), ASPET-Astellas Award, a PhRMA Foundation “Center of Excellence in Clinical Pharmacology” Award, and the Gerstner Family Career Development Awards in Individualized Medicine.
The authors declare that they have no competing interests.
LL, NN and LW designed the study and wrote the manuscript. LL and NN performed the experiments. LL, BF, KK, GJ, AB and LW analyzed the data. LL, NN, BK and LW wrote the manuscript. All authors read and approved the final manuscript.
Electronic supplementary material
Additional file 1:Contains supplementary figures and tables. Figure S1. Imputation analysis for top SNP loci associated with the response of (A) gemcitabine, (B) AraC, and (C) both drugs. Figure S2. Validation of imputed genotypes. The x-axis indicates actual genotype by TaqMan assay. The y-axis represents imputed genotype, which was estimated as the count of a particular allele. The squared difference between the imputed genotype and actual genotype was calculated based on counting the same allele. Avg sq difference = average squared difference. Figure S3. Effect of FKBP5 on cellular sensitivity to either gemcitabine or AraC in which PIGB was stably knock down cell lines. Table S1. Top SNPs that were associated with (A) Gemcitabine, (B) AraC, or (C) both drugs with P values <10-3 during the GWAS. Table S2. SNPs associated with both expression and cytotoxicity data for either (A) Gemcitabine or (B) AraC from the “integrated analyses” with P values <10-3. Table S3. The top 9 loci for Gemcitabine (A) and the top 9 loci for AraC (B) that were associated with drug response-IC50 values using our previous 550,000 SNP array data. (PDF 437 KB)
Authors’ original submitted files for images
Below are the links to the authors’ original submitted files for images.
About this article
Cite this article
Li, L., Fridley, B.L., Kalari, K. et al. Discovery of genetic biomarkers contributing to variation in drug response of cytidine analogues using human lymphoblastoid cell lines. BMC Genomics 15, 93 (2014). https://doi.org/10.1186/1471-2164-15-93
- Cytidine analogues
- Cytosine arabinoside
- Lymphoblastoid cell lines
- Expression array
- Genome-wide SNPs
- Genome-wide association study
- Functional genomics
- Translational research