Profile of the breast cancer susceptibility marker rs4245739 identifies a role for miRNAs.

OBJECTIVE: To determine the influence of the single nucleotide polymorphism (SNP) rs4245739 on the binding and expression of microRNAs and subsequent MDM4 expression and the correlation of these factors with clinical determinants of ER-negative breast cancers.
METHODS: FindTar and miRanda were used to detect the manner in which potential microRNAs are affected by the SNP rs4245739-flanking sequence. RNA sequencing data for ER-negative breast cancer from The Cancer Genome Atlas (TCGA) were used to compare the expression of miR-184, miR-191, miR-193a, miR-378, and MDM4 in different rs4245739 genotypes.
RESULTS: Comparison of ER-negative cancer patients with and without the expression of miR-191 as well as profile microRNAs (miR-184, miR-191, miR-193a and miR-378 altogether) can differentiate the expression of MDM4 among different rs4245739 genotypes. Although simple genotyping alone did not reveal significant clinical relationships, the combination of genotyping and microRNA profiles was able to significantly differentiate individuals with larger tumor size and lower number of involved lymph nodes (P < 0.05) in the risk group (A allele).
CONCLUSIONS: We present two novel methods to analyze SNPs within 3'UTRs that use: (i) a single miRNA marker expression and (ii) an expression profile of miRNAs predicted to bind to the SNP region. We demonstrate that the application of these two methods, in particular the miRNA profile approach, permits detection of new molecular and clinical features related to the rs4245739 variant in ER-negative breast cancer.

Introduction

Genetic variations, including single nucleotide polymorphisms (SNPs), have emerged as an important resource in the understanding of carcinogenesis[1-3]. Recently, several studies have reported the role of SNPs in carcinogenesis through their associations with the binding and expression of microRNA (miRNAs)[4-6]. MiRNAs are short non-coding RNAs (~18–25 nucleotides) that negatively regulate gene expression post-transcriptionally by directly binding to mRNA targets[7,8]. MiRNAs have been associated with the pathogenesis of almost every human cancer through regulation of several important biological pathways, including cell proliferation, migration, and apoptosis[8-11]. Unique patterns of miRNA expression have been proposed as clinical biomarkers and potential therapeutic targets in several cancers[12-15].

One SNP, rs4245739, has been identified as a risk factor for several cancers, including esophageal cancer[16], non-Hodgkin’s lymphoma[17], prostate cancer[18], lung cancer[19], and ovarian cancer[20]. Rs4245739(A > C) has also been associated with susceptibility to breast cancer [21], the most common cancer among women and the fifth leading cause of cancer mortality worldwide[22,23]. In a study involving a total of more than 10,000 breast cancers negative for estrogen receptor (ER) expression and 75,000 control cases from patients of Caucasian ancestry, pooled results showed that rs4245739 was associated with the risk of ER-negative breast cancer (P < 0.01) but not the risk of ER-positive breast cancer, with an odds ratio (OR) for the A risk allele of 1.14 (95% CI = 1.10–1.18) [22]. Subsequent studies from China (1,100 breast cancer cases/1,400 controls) and Norway (1,717 breast cancer patients/1,870 controls) have also documented rs4245739 as a risk variant for breast cancer[21,24].

Rs4245739(A > C) is located in the three prime untranslated region (3’UTR) of MDM4[25], an oncogene that negatively regulates p53 expression by directly binding and masking the transactivation domain in response to DNA damage[26]. This SNP creates a target binding site for the miRNA miR-191-5p, which has been linked to several cancer types[18,25]. A functional study revealed miRNA-mediated fine-tuning of MDM4 expression as a result of miR-191-5p binding to the C-allele but not the A-allele of rs4245739[25]. Therefore, in A-allele genotypes, MDM4 is overexpressed and associated with a higher risk of breast cancer[25]. In ovarian cancer, patients with the AA genotype showed significantly increased risk of recurrence and mortality[25]. Upregulation of miR-191-5p has been reported in breast cancer, although its specific role in carcinogenesis or progression in this cancer type is not yet known[27].

In this study, we used a combination of genetic variant data, in silico prediction, and miRNA expression profiles to demonstrate processes by which molecular mechanisms involved in the SNP located at the 3’UTR of the MDM4 gene can alter the expression of several miRNAs associated with specific clinical features, i.e., tumor size and lymph node infiltration in breast cancer.

Materials and methods

Tumor array and sequencing data

Affymetrix Genome-wide SNP 6.0 microarray (Affymetrix Inc., Santa Clara, CA, USA), gene-level miRNA, and gene-level mRNA sequencing data were retrieved for The Cancer Genome Atlas (TCGA) breast cancer cohort[28]. RNA sequencing data determined with RNA-Seq by Expectation-Maximization (RSEM) software were used. ER-negativity was determined from immunohistochemical records. Altogether, SNP array and expression data were available for 189 primary ER-negative breast cancers. For these cancers, we were able to find corresponding clinical stage, tumor size, and lymph node status.

SNP survey

SNPs listed on dbSNP were extracted, and those with a minor allele frequency (MAF) of less than 0.1 were removed to yield a total of 11,052,394 SNPs[29]. The genomic location of SNPs was determined with reference to the human genome assembly GRCh37.

SNP genotyping

Genotypes for the 906,600 SNP probes of the Affymetrix SNP 6.0 platform (Affymetrix, Inc.) were determined by first mapping probe set IDs to dbSNP rs IDs. The R package genomewidesnp6Crlmm was then used for genotype calls, with calls having a confidence level of less than 90% filtered out. All 189 samples were retained for rs4245739.

MiRNA-mRNA binding alteration prediction

We used the algorithms FindTar (http://bio.sz.tsinghua.edu.cn/findtar)[30] and miRanda (http://www.microrna.org/microrna/home.do) to predict the mRNA targets of miRNAs.

Bioinformatic and statistical analysis

Bioinformatic and statistical analyses were conducted with the R statistical package version 3.1.2. Cut-offs for the different methods were determined as described previously[31]. P values less than 0.05 were considered statistically significant unless otherwise stated.

Results

Screen of SNPs located at UTRs and miRNA-mRNA binding alterations due to rs4245739

We analyzed the role of SNPs present in the 3’UTRs of the human genome. After discarding those with a MAF of less than 0.1, we obtained 11,052,394 SNPs, the majority of which were inter-genic. Of the SNPs located in genes, more than 100,000 were in 3’UTRs, a quantity almost two and four times the amounts in the coding regions (non-synonymous only) and 5’UTRs, respectively. We next calculated the probability of finding SNPs in these regions calibrated by their total base length in the genome (, left). This probability was higher for 3’UTR regions compared to the coding regions and 5’UTRs, suggesting a selective pressure for this occurrence. We therefore predicted that some of the 3’UTR-located SNPs may alter the binding ability of regulatory miRNAs (, right), which may in turn mediate allele-dependent mRNA post-transcriptional stability. We assessed this phenomenon by analyzing rs4245739, which is present in the 3’UTR of the MDM4 gene and for which the A allele predisposes individuals to ER-negative breast cancer ().[22] We used the FindTar3[30] and miRanda[32] algorithms to predict miRNA-mRNA binding alterations between risk and wild type alleles of 3’UTRs produced by rs4245739 (), which showed not only miR-191 but other miRNAs to be affected by this SNP.

MiRNA binding alterations due to SNP rs4245739. (A) Distribution of SNPs located at the 3’UTR and 5’UTR per 100,000 bases in comparison to coding regions (NS) in the human genome (left panel). Schematic description of how SNP located at the 3’UTR might influence binding sites of miRNA (right panel). (B) Schematic representation of SNP rs4245739, which is located at the 3’UTR of the MDM4 gene. (C) Potential binding alterations of several miRNAs due to rs4245739 that affect not only miR-191 but also 16 additional miRNAs.

Expression of miRNAs affected by SNP rs4245739 genotypes in patients with ER-negative breast cancer

We used an ER-negative breast cancer data set retrieved from the TCGA to study the manner in which the potential alterations in miRNA binding result in phenotypic changes. First, we obtained the rs4245739 genotype of individuals from this data set and observed the homozygotic genotype (CC) and heterozygotic genotype (AC) in approximately 20% of primary tumors (). Interestingly, there were very few genotypic changes between germline and tumor tissue, with only conversion from the heterozygote to either the homozygote C-allele or homozygote A-allele being observed. We next focused on the expression of miRNAs with binding to MDM4 that could be altered in an allele-dependent manner (). Only miR-184, miR-191, miR-193a, and miR-378 had an appreciable average expression, albeit with a large variance. Notably, two of these miRNAs (miR-193a and miR-378) were predicted to bind to the A-allele MDM4 3’UTR only, whereas the other two (miR-184 and miR-191) were predicted to bind to it only in the presence of the C-allele. Examination of whether the expression of these miRNAs differed between genotypes did not yield significant results ().

Expression of miRNAs according to different rs4245739 genotypes in ER-negative breast cancer. (A) Distribution of ER-negative breast cancer patients according to rs4245739 genotypes. (B) Expression of 17 miRNAs with binding to MDM4 mRNAs that is altered owing to rs4245739. (C) Expression of miR-184, miR-191, miR-193a, and miR-378 in different rs4245739 genotypes.

Characterization of rs4245739-affected miRNA expression and 3’UTR SNP functions

Reasoning that the expression level of these miRNAs may complement our ability to understand 3’UTR SNP function and its effect on disease, we approached the characterization of rs4245739 in three different ways (). In the first approach, we simply grouped tumors according to rs4245739 genotype (, left). In the second method, we additionally stratified patients according to expression or lack of expression of miRNA-191 (, center), which would enable us to assert with more confidence whether allele-dependent miRNA modulation of MDM4 was occurring. Our final approach for classifying samples consisted of generating an miRNA profile of each individual (, right). In the case of rs4245739, we retrieved from the TCGA data sets the expression levels of the four miRNAs predicted to affect MDM4 expression (miR-184, miR-191, miR-193a, and miR-378). Given that miR-193a and miR-378 would bind to the A-allele-containing 3’UTR, and miR-184 and miR-191 would bind to the C-allele-containing 3’UTR, we divided the individuals into an A-allele SNP+/miR-193a-/miR-378- group and a C-allele SNP+/miR-184+/miR-191+ group. These two groups would represent A-genotype individuals without miRNA-mediated MDM4 alterations, and C-genotype individuals with miRNA-mediated MDM4 alterations. We compared these two novel approaches to 3’UTR SNP analysis with the genotype-only classification by performing differential MDM4 gene expression analysis between the groups with each approach. In order to categorize tumors into miR-191-expressing and non-miR-191-expressing cases, we compared the results obtained with different thresholds against the results from the simple genotyping approach (, left). We did not observe significant changes in MDM4 expression among any of the groups generated by simple genotyping or the miR-191-expressing approach. Nevertheless, when we applied the miRNA-profile approach with different miRNA expression thresholds (differences in , right, and P values in , left), we observed significant changes in MDM4 expression. Unexpectedly, we found overexpression of MDM4 in the case of miRNA binding to its 3’UTR. In light of these results, we analyzed the expression of other genes located near (< 2Mb) the rs4245739 SNP using the miRNA profile approach. Three other genes (MFSD4, SLC26A9, and SLC41A1) were significantly repressed in the miRNA-bound groups (, right).

Characterization of rs4245739-affected miRNA expression using several methods to show 3’UTR function.

Activated biological pathways in patients with ER-negative cancer with different rs4245739 genotypes and miRNA expression profiles

We next analyzed the full transcriptome to detect whether gene expression differences were conserved among the different methods. We observed 100% overlap between the sets of differentially expressed genes identified by these two techniques up to a threshold of 300 units of miR-191 expression, at which point the overlap rapidly fell below 50% (, left). At the same time, the number of genes detected by the marker method remained relatively stable across different expression thresholds (, right). Finally, when we analyzed the results obtained with the third method, which combined the expression levels of miRNA predicted to bind to a given allele to give an integrated score, we found that the differentially expressed gene overlap between the combined approach and simple genotyping approaches rapidly decreased (, left). As with the previous method comparison, differences were not related to the number of genes obtained with each approach (, right). To analyze the biological importance of the genes obtained, we used the 95th percentile of the entire miRNA expression data set in the case of the miR-191 marker. As a result, we obtained two groups: one comprising individuals with the C-genotype and miR-191 expression above the aforementioned cut-off (n = 37), and another comprising all other individuals (n = 157). In the case of the miRNA profile method, we restricted our samples to the 70% of individuals following the miRNA profile results. The two groups comprised 75 individuals with the A-genotype, and 45 individuals with the C-allele. We performed pathway analysis of these results using the DAVID webserver (). We found processes such as “immune response” and “cell adhesion” to be overrepresented in all three methods, but each method also detected specific pathways. For example, a protocadherin cluster was found only when using the miRNA profile.

Clinical relevance of characterization of SNP rs4245739 genotypes and miRNA profiles

Finally, we assessed the strength of association between the results of these new methods and clinical parameters. The use of simple genotyping alone failed to reveal any significant relationships ( and and , P values marked by dotted lines). When we examined tumor size categories of samples, we found significant differences among groups using both methods at higher thresholds (P values in , differences in ). In addition, the number of lymph nodes affected in individuals was found to be significantly reduced in the C-allele group obtained by the miRNA profile method (P values in , right; differences in Supplementary , right), using several different thresholds. Results with the marker method, similar to those using the simple genotyping approach, did not reach significance in this association (, left; Supplementary , left). Associations with other clinical parameters, including breast cancer mortality (Supplementary ) and clinical stage (Supplementary ), were not statistically significant. Our novel methods thus improve the prediction of these associations compared to mere use of the presence or absence of a SNP variant.

Discussion

We demonstrated the use of a single miRNA expression marker and expression profile of miRNAs associated with rs4245739, a SNP located at the 3’UTR of the MDM4 gene, to identify tumor features linked to this variant. Our methods, in particular, use of the expression profile, were superior compared to classification based on a single SNP in predicting disease characteristics and outcomes among patients with ER-negative breast cancer. In addition, we found overexpression of MDM4 in the case of miRNA binding to its 3’UTR. Previous research has shown that rs4245739 affects binding of miR-191 and regulates expression of MDM4 mRNA in ovarian cancer[25]. These findings demonstrate the complexity of the processes governing miRNA-mRNA interactions[18,33]. We combined SNP genotyping data with miRNA as well as MDM4 expression in ER-negative breast cancer. Single miRNA expression in the presence of the SNP variant did not influence MDM4 expression (, left). However, using the miRNA profile approach, we were able to observe differential changes in MDM4 expression (, right; , left). Our results confirm the findings of a previous study of ovarian cancer by Wynandaele and colleagues[25].

More importantly, by combining SNP variant presence and the miRNA profile approach, we were able to differentiate subgroups of patients at risk with greater tumor size and number of lymph nodes infiltrated, indicating the potential for the use of SNPs in this way as diagnostic and prognostic markers in ER-negative breast cancer. Tumor size is used together with nodal and metastasis status to determine clinical breast cancer staging and suitable breast cancer management. In addition, around 30% of breast cancers are ER-negative tumors, with a higher proportion in young breast cancer patients, and these tumors constitute a subgroup with relatively poorer prognosis and fewer treatment options[34-36]. Our enrichment results showed a common theme of immune response across the three different methods; these findings that are consistent with previous results highlighting the importance of immune pathways in ER-negative breast cancer[37-39].

We found differential MDM4 expression using microRNA profiles in relation to different genotypes. Our findings support the previous notion that SNPs located at 3’UTRs may lead to dysregulation of target mRNAs and are potentially implicated in the development of cancer and other chronic diseases[25]. Indeed, many studies have shown that SNP variants in the 3’UTR elevate the risk of developing several cancer types. Saetrom and associates[40] showed that rs1434536 located at the bone morphogenetic receptor type 1B gene (BMPR1B) 3’UTR influences the binding of miR-125b and correlates with breast cancer risk. Located at the estrogen receptor-1 gene (ESR1) 3’UTR, rs9341070 is a risk variant for breast cancer. The SNP affects interaction with miR-206 and influences the expression of ESR1[41]. As a risk variant for squamous cell carcinoma of the head and neck, SNP rs8126 (T > C) is located at the TNFAIP2 3’UTR and influences the expression of TNFAIP2[42]. In addition, rs3134615, which is a susceptibility variant for small-cell lung cancer, is located at the MYCL1 3’UTR and influences the affinity for miR-1827[43]. Rs34764978, which is mapped at the 3’UTR of human dihydrofolate reductase (DHFR), affects the binding and expression of miR-24 and correlates significantly with resistance to the chemotherapeutic agent methotrexate[44]. These data show that SNPs of specific genes located at the 3’UTR play an important role in the pathogenesis of several diseases through modulating the binding of particular miRNAs.

MiR-191[25] as well as miR-184, miR-193, and miR378 interact and modulate MDM4 expression to certain rs4245739 genotypes. MDM4 is an oncogene that inhibits p53 tumor suppressor proteins during regulation of cell proliferation and response to DNA damage[25]. Nagpal et al[27]. have reported that miR-191 is an estrogen-sensitive microRNA that acts in regulation of cell proliferation through inhibition of BDNF, CDK6, and SATB1. MiR-184 regulates the PI3K/AKT/mTOR pathway during mammary development and tumor progression[45] as well as p53/p21 activity during the cell cycle and cell migration[46]. MiR-193a inhibits cell growth and migration in breast cancer cells by directly targeting NLN, CCND1, PLAU, and SEPN1[47]. MiR-378 induces a metabolic shift from the glycolytic to the oxidative pathway through PGC-1β/ERRγ[48], and modulates therapeutic responses to tamoxifen by targeting GOLT1A[49]. In addition, gene enrichment pathway analysis revealed that biological processes such as “immune response” and “cell adhesion” were overrepresented in all three methods we used, suggesting that these microRNAs play roles in immune response and cell migration. As its binding affinity for the MDM4 3’UTR is affected by SNP rs4245739, miR-92b is an important candidate for further functional studies because the precursor sequence forms a stem loop for miR-17-92 cluster, which plays a significant role in the carcinogenesis of many cancers, including breast cancer[50]. MiR-17-92 cluster is known as oncomir-1[50,51]. Further investigations could involve experiments to confirm functional rs4245739, miR-92b, and MDM4 or TP53 expression with the use of luciferase reporter assays and overexpression of miR-92 mimics in combination with reverse transcription-PCR in cell lines.

Acknowledgements

Wahyu Wulaningsih was supported by grants from UK Medical Research Councils (Grant No. MC_UU_12019/2 and MC_UU_12019/4). Sumadi Lukman Anwar received a grant from PTUPT (Grant No. Ristekdikti 09_18).

Conflict of interest statement

No potential conflicts of interest are disclosed.

1

DAVID results for rs4245739 analysis using method 1 (only SNP presence/absence)

Term	Count	Fold Enrichment	P
IL12/Stat4 signaling pathway in TH1	5	29.9	0.0001
Signal peptide	38	1.8	0.0002
Cell activation	9	4.6	0.0007
Cell adhesion	14	2.9	0.0008
ER	68	1.3	0.0031
Immune response	13	2.7	0.0024

2

DAVID results for rs4245739 analysis using method 2 (SNP and miR-191 expression)

Term	Count	Fold Enrichment	P
Immune response	114	1.7	< 0.0001
Cell adhesion	24	4.6	< 0.0001
Extracellular matrix	20	4.2	< 0.0001
Glycoprotein	86	1.6	< 0.0001
EGF-like domain	13	4.5	0.0001
Focal adhesion	10	4.4	0.0003
Cell morphogenesis	14	2.9	0.0012
Negative regulation of proliferation	13	2.6	0.0041
Inhibition of matrix metalloproteinases	3	26.6	0.0045

3

DAVID results for rs4245739 analysis using method 3 (SNP and miR profile)

Term	Count	Fold Enrichment	P
Protocadherin beta	5	33.6	< 0.0001
Basolateral plasma membrane	11	6	< 0.0001
Cell adhesion	15	4.1	< 0.0001
Immune response	23	2.6	< 0.0001
Spectrin/alpha-actinin	5	18.4	0.0001
Integrin binding	6	10.7	0.0002
ARVC (cardiomyopathy)	5	8.4	0.0026
Epidermolysis bullosa	3	29.0	0.0046