Introduction
Breast cancer (BC) is the most prevalent type of cancer affecting women, and the second most common type of cancer as a whole. Globally, the incidence of BC is one in nine women (1). Despite substantial advancements in diagnosis and treatment, the mortality rates for BC are still at 15% (2).
One of the most aggressive variants of BC is triple-negative BC (TNBC). This subgroup of BC cells does not express the receptors for estrogen (ER), progesterone receptor (PR) and human epidermal growth factor 2 (HER2). The absence of these proteins renders commonly used clinical interventions, such as inhibiting aromatase and blocking hormone receptors ineffectual for TNBC therapy (3,4). Treatment is therefore limited to surgical resection and adjuvant chemotherapy, making TNBC a significant unmet clinical need (5).
Synthetic glucocorticoids, such as dexamethasone are commonly used in tumor therapy (6). However, recent data have indicated that mortality is increased in patients with TNBC variants overexpressing glucocorticoid receptor (GR). GR signaling activates oncogenes, inhibits apoptosis and represses tumor suppressor genes in TNBC, leading to unfavorable clinical outcomes (7). A high GR expression is also associated with therapy resistance and increased recurrence (8).
In addition to regulating protein-coding genes, GR is known to activate non-coding RNAs, including microRNAs (miRNAs or miRs) (9). These molecules of 18-22 nt in length modulate cellular gene expression by specifically binding complementary mRNA sequences and repressing their translation. Highly malignant tumor cells commonly display dysregulated miRNA profiles, leading to oncogenic and anti-apoptotic signaling (10). For TNBC, previous studies have discovered signatures of altered miRNA expression that distinguish cancer cells from surrounding tissues, and predict the receptor status for ER, PR and HER2 (11,12). Additionally, the analysis of miRNAs highlights the heterogeneity of TNBC phenotypes and their respective signaling pathways (13). Unraveling miRNA regulation for individual TNBC variants is therefore crucial to understand pathogenesis and draw consequences about potential therapeutic interventions.
Recent studies have indicated that TNBC modulates distant cells by secreting signaling factors, including extracellular vesicles (EVs). EVs shed from tumor cells are enriched in specific miRNAs that might contribute to tumor progression and metastasis (14). Indeed, BC EVs have been found to induce the proliferation, migration and invasion of recipient cells, thus enhancing disease progression (15–17). Exploiting circulating miRNAs for the diagnosis and prognosis of TNBC patients might yield powerful biomarkers that are easily accessible via liquid biopsy (18).
To the best of our knowledge, however, to date, there is no study available investigating the effects of GR signaling on miRNA regulation in TNBC. Thus, in this study, to address this knowledge gap, we analyzed cellular and vesicular miRNA profiles in an in vitro model of GR-overexpressing TNBC by high-throughput next-generation sequencing (NGS). Surprisingly, no statistically significant alterations in EV miRNAs were detected in two TNBC cell lines upon GR overexpression. We did, however, detect a small set of cellular miRNAs regulated by GR in a cell line-specific manner. As the validation of miRNA target genes yielded ambiguous results, we concluded that the unfavorable influence of a high GR expression on TNBC phenotypes is not mediated by miRNAs to a significant extent.
Materials and methods
TNBC cell culture and transfection
The human TNBC cell lines, MDA-MB-231, MDA-MB-436 and MDA-MB-468, were purchased from the Leibnitz Institute DSMZ-German Collection of Microorganisms and Cell Culture (Braunschweig, Germany) and Cell line services (Eppelheim, Germany).
The cells were cultured in T75 flasks in a monolayer in DMEM (Sigma-Aldrich, Hamburg, Germany) containing 4.5 g/l glucose, 1% L-glutamine, 10% exosome-depleted fetal bovine serum (FBS) (BioCat GmbH, Heidelberg, Germany), 160 ng/l cortisol and 1% penicillin/streptomycin (Invitrogen, Karlsruhe, Germany). Cultures were maintained in a humidified atmosphere at 37°C and 5% CO2. For the experiments studying cellular miRNAs, the MDA-MB-231, MDA-MB-436 and MDA-MB-468 cells were seeded 4.5 h prior to transfection at a concentration of 1×105 cells per well in 24-well plates (Greiner Bio-One, Frickenhausen, Germany) using 0.5 ml culture medium.
GR overexpression was induced by transfecting the TNBC cells with nuclear receptor subfamily 3 group C member 1 (NR3C1)-encoding DNA plasmids. The coding sequence (CCDS4278.1) of the predominant glucocorticoid receptor gene NR3C1 transcript variant 1 (NM_000176.2) was synthesized with the restriction sites KpnI and XhoI on the 5′ and 3′ end, respectively (MWG Eurofins, Ebersberg, Germany), and cloned in frame into the pcDNA6/V5-His A vector (Life Technologies, Darmstadt, Germany). The cells were transfected with the NR3C1 plasmid (0.4 µg) for 24 h using Lipofectamine 2000 transfection reagent (Invitrogen). Following 30 h of cultivation, the cells had reached 90% confluency, and were rinsed once with HBSS before proceeding to total RNA extraction.
For experiments studying EV miRNAs, the MDA-MB-231 and MDA-MB-468 cells were seeded at a concentration of 3×106 cells per well in 6-well plates (Greiner Bio-One). Both parental and transfected cells were seeded in 3 wells with 2.5 ml culturing medium each. Transfection was performed with 2.5 µg of plasmid harboring the coding sequence of NR3C1.
For all experiments, untransfected cells with endogenous GR expression were used as control samples. Three independent technical replicates per cell line were analyzed.
Validation of NR3C1 overexpression
To evaluate the effectiveness of transfection, the NR3C1 mRNA levels were quantified by reverse transcription-quantitative (real-time) PCR (RT-qPCR). Total RNA from each cell line was initially reverse transcribed using the QuantiTect Reverse Transcription kit (Qiagen, Hilden, Germany) according to the manufacturer's instructions. Quantitative PCR was performed using the Prime PCR Assay NR3C1 and GAPDH human SsoAdvanced universal supermix (Bio-Rad, Munich, Germany) and 10 ng of template cDNA. PCR reactions were run on a MiniOpticon real-time PCR system (Bio-Rad). Additionally, the transcription levels of NR3C1 downstream targets dual specificity phosphatase 1 (DUSP1; NM_004417), serum/glucocorticoid regulated kinase 1 (SGK1; NM_005627.3) and glucocorticoid-induced leucine zipper protein (GILZ; NM_198057.2) were assessed by RT-qPCR. For each cell line, total RNA from three biological replicates was reverse transcribed using the QuantiTect Reverse Transcription kit (Qiagen). Subsequently, 8 ng cDNA were analyzed in a 10 µl reaction volume of SsoFast EvaGreen Supermix (Bio-Rad) and 300 nM primers (Sigma-Aldrich). Real-time PCR was carried out on triplicate samples using a Rotor-Gene Q thermal cycler (Qiagen) and a thermal profile for polymerase activation (95°C for 1 min) and 45 cycles of amplification (95°C for 10 sec 60°C for 15 sec, 65°C for 45 sec). The primer sequences are provided in Table I. The expression of NR3C1 and its downstream targets was normalized to GAPDH, a stable reference gene for breast cancer cells (19). Relative quantification was carried out using the ΔΔCq method (20). Statistical significance was determined using the Student's t-test. Values of P<0.05 were considered to indicate statistically significant differences.
Extracellular vesicle isolation and characterization
For EV isolation, 7.5 ml of cell culture supernatant were collected from the parental and transfected cells after 30 h of cultivation, and centrifuged (3,200 × g, 5 min) to remove the cellular debris. EVs were isolated from pre-cleared supernatant using the miRCURY Exosome Isolation kit - Cells, urine and CSF according to the manufacturer's instructions (Exiqon, Vedbaek, Denmark). EV pellets were resuspended in either lysis buffer for RNA extraction, or PBS for vesicle characterization.
For nanoparticle tracking analysis (NTA), the EVs were diluted in particle-free PBS and analyzed on a NanoSight LM10 (Malvern Instruments GmbH, Herrenberg, Germany) using a 408 nm laser and NTA 3.0 software. Four videos of 30 sec each were captured, and analyzed using default settings for blur and minimum track length, and a detection threshold of two.
RNA extraction and NGS library preparation
Total RNA was isolated from the cells and EVs using the miRCURY RNA Isolation kit – Cell and Plant (Exiqon,) according to the manufacturer's instructoins. Cellular RNA was quantified using a nanophotometer (Implen GmbH, Munich, Germany), and RNA integrity was assessed by capillary electrophoresis on the Bioanalyzer 2100 using the RNA 6000 Nano kit (Agilent Technologies, Waldbronn, Germany). EV RNA was analyzed using the Agilent Small RNA kit (Agilent Technologies).
Sequencing libraries were constructed from 190 ng of cellular RNA, or the entire EV RNA isolated from 7.5 ml conditioned culture medium, respectively. Library preparation was performed as previously described by Spornraft et al (21). Briefly, the RNA was adaptor-ligated, reverse-transcribed, amplified by PCR and barcoded using the NEBNext Multiplex Small RNA Library Prep Set for Illumina (New England BioLabs Inc., Frankfurt, Germany). Adaptors and primers were diluted 1:2 in nuclease-free water to accommodate the low RNA input. Size selection of pooled PCR products was performed by agarose gel electrophoresis (4%), cutting out bands with 130 to 150 bp fragments. The purity and concentration of the libraries extracted from the gel were verified by capillary electrophoresis using the High Sensitivity DNA kit on the Bioanalyzer 2100 (Agilent Technologies). Finally, the libraries were subjected to Illumina single-end sequencing-by-synthesis using 50 cycles on the HiSeq 2500 (Illumina Inc., San Diego, CA, USA).
Data processing and differential gene expression analysis
FastQC (version 0.11.5) was used to assess the sequence length distribution and quality of the NGS data, as previously described (22). Adaptor sequences were trimmed using BTRIM, and all reads without adaptors were discarded (23). Additionally, reads shorter than 15 nt were excluded from the data set (24). Prior to miRNA analysis, reads pertaining to ribosomal RNA (rRNA), transfer RNA (tRNA), small nuclear RNA (snRNA) and small nucleolar RNA (snoRNA) were removed by mapping to sequences obtained from RNAcentral (25). The remaining reads were then aligned to miRBase (version 21) (26). Mapping was carried out using Bowtie and the 'best' alignment algorithm, allowing one mismatch for both RNAcentral and miRBase (27). Final read count tables were generated by sorting and indexing aligned reads using SAMtools, and calling the sum of hits per miRNA sequence (28). Differential gene expression analysis was subsequently performed via the bioconductor package DESeq2 (version 1.8.1), using the Benjamini-Hochberg method to correct for false discovery (29). A log2 fold change ≥|1| and an adjusted p-value (Padj) of ≤0.05 were set as thresholds to identify significantly regulated miRNAs. Only miRNAs with a mean expression of at least 50 counts were included in the analysis. Principal component analysis (regularized log-transformed, sizefactor-corrected counts obtained from DESeq2), and data visualization were performed in R (version 3.4.0) using the packages gplots, ggfortify, genefilter and RColorBrewer.
Validation of regulated miRNAs
Based on the NGS data, differentially regulated cellular miRNAs were validated by RT-qPCR. First, 111 ng of RNA were reverse transcribed in triplicate using the miScript II RT kit (Qiagen) according to the manufacturer's instructions. A total of 1 µl cDNA was subjected to real-time PCR in a 10 µl reaction volume using the miScript SYBR-Green PCR kit (Qiagen). Reactions were run on a CFX384 real-time PCR detection system (Bio-Rad) using the recommended protocol of polymerase activation (95°C for 15 min) and 45 cycles of amplification (94°C for 15 sec, 55°C for 30 sec, 70°C for 30 sec). Quantification cycle (Cq) values were determined automatically using default threshold settings, and Cq values above 37 were manually set to 40. The NGS data was utilized to assess potential reference miRNAs using the geNorm and NormFinder algorithms (30,31). Cq values of regulated miRNAs were subsequently normalized to the geometric mean of the following reference miRNAs: miR-24-3p, miR-25-3p and miR-148b-3p for MDA-MB-231 and MDA-MB-436, and miR-24-3p, miR-25-3p and let-7a-5p for MDA-MB-468. Relative quantification was carried out using the ΔΔCq method (20). Statistical significance was assessed using Student's t-test. Values of P<0.05 were considered to indicate statistically significant differences.
Prediction and quantification of miRNA target genes
miRWalk 2.0 was used to predict mRNAs targeted by miR-203a-3p (32). Four target genes known to be associated with metastasis in solid tumors [Actin, gamma 2, smooth muscle, enteric (ACTG2), calponin 1 (CNN1), major histocompatibility complex, class II, DP beta 1 (HLA-DPB1) and myosin light chain kinase (MYLK)] were selected for analysis by RT-qPCR (33). The expression of these genes was quantified in the same cellular MDA-MB-436 samples previously used for NGS. Initially, RNA was reverse transcribed in triplicate using the QuantiTect Reverse Transcription kit (Qiagen) according to the manufacturer's instructions. Quantitative PCR was then performed using Prime PCR Assays, SsoAdvanced Universal Supermix (Bio-Rad) and 10 ng of template cDNA. PCR reactions were run on a MiniOpticon real-time PCR system (Bio-Rad). Target gene Cq values were normalized to GAPDH (19).
Results
Plasmid transfection induces the expression of NR3C1 and downstream effectors
The transfection of the TNBC cells with NR3C1-coding plasmids was validated by RT-qPCR (Fig. 1A). GR expression was significantly increased in the transfected MDA-MB-231 (111-fold, P=1.75E-4), MDA-MB-436 (172-fold, P=2.51E-4) and MDA-MB-468 (335-fold, P=6.39E-6) cells.
To assess overexpression-induced changes in GR signaling, we additionally quantified the expression levels of GR target genes DUSP1 (Fig. 1B), SGK1 (Fig. 1C) and GILZ (Fig. 1D). While DUSP1 and SGK1 were significantly upregulated in all cell lines, GILZ expression was not altered significantly in the MDA-MB-468 cells (P=0.51). In the MDA-MB-231 and MDA-MB-436 cells, however, GILZ expression was significantly increased (P=2.79E-3 and P=3.79E-7, respectively).
TNBC secretes extracellular vesicles carrying RNA
EVs isolated from the conditioned media of the parental MDA-MB-231 and MDA-MB-468 cells were characterized by NTA (Fig. 2). Single-particle analysis revealed a narrow size distribution with mean particle diameters of 119.0±74.4 nm (mode, 87.5 nm) and 140.0±116.7 nm (mode, 117.5 nm) for the MDA-MB-231 and MDA-MB-468 cells, respectively. Despite the similarities in diameter, significantly more particles were isolated from the MDA-MB-231 cells (P<0.001).
Prior to sequencing, RNA extracted from EV preparations was analyzed by capillary electrophoresis. Samples from both cell lines were found to be enriched in small RNA species <150 nt without obvious differences in size profiles between the parental and transfected cells. Full electropherograms for small RNA analysis are provided in Fig. 3.
GR overexpression does not alter vesicular miRNA profiles
EV RNA from the MDA-MB-231 and MDA-MB-468 cells was profiled by small RNA-Seq. The mean per-replicate library sizes are provided in Table II. After mapping to miRBase, between 736 (MDA-MB-468, transfected) and 796 (MDA-MB-231, parental) distinct miRNA transcripts were detected in at least one sample. For EVs from both cell lines, there was significant overlap in the 10 most highly expressed miRNAs between the treatment groups (Table III).
Differential expression of miRNAs was assessed in EVs from the parental and transfected MDA-MB-231 and MDA-MB-468 cells. While individual cell lines were clearly distinguished by principal component analysis (Fig. 4), the overexpression of GR did not lead to noticeable changes in miRNA expression. None of the miRNAs detected in small RNA-Seq displayed statistically significant regulation between endogenous and artificially induced GR expression.
Changes in cellular RNA profiles upon GR overexpression
Cellular RNA was initially analyzed by capillary electrophoresis to assess its suitability for NGS analysis. For all cell lines, samples from both parental and transfected cells featured excellent RNA integrity, as indicated by the RNA integrity number (RIN) values >9. Bioanalyzer electropherograms for cellular RNA are shown in Fig. 5.
In the NGS data, both the mean size of sequencing libraries and the number of detected miRNAs were higher than in the EV samples (Table IV compared with Table II). The most highly expressed miRNAs in all of the three parental cell lines displayed a high degree of similarity, sharing 8 of the top 10 miRNAs (Table V). Similarly, 7 of the top 10 most highly expressed miRNAs were common to all transfected cell lines.
Differential gene expression analysis revealed slight, yet statistically significant changes in specific miRNAs during GR overexpression (Table VI). Of note, a different set of GR-responsive miRNAs was detected in each of the TNBC cell lines studied herein, highlighting the heterogeneity of molecular signaling. As shown in Fig. 6, miRNA expression clearly separated individual cell lines, but exhibited only minor differences between the transfected and parental cells. Based on the expression profiles of all miRNAs, principal component 3 (PC3), distinguished parental and GR-overexpressing cells (Fig. 7A). Limiting the input for analysis to the 500 miRNAs with highest variance, however, a reduced separation of groups was observed (Fig. 7B).
We then assessed differential miRNA regulation between endogenous and induced GR expression using RT-qPCR and the Student's t-test. Of the 7 miRNAs found to be significantly regulated in the NGS data, only miR-203a-3p was validated with statistical significance. In the transfected MDA-MB-436 cells, it was upregulated with a log2 fold change of 0.63 (Fig. 8).
To assess the potential biological functions of miR-203a-3p in the MDA-MB-436 cells, in which it was found to be significantly increased in, we quantified the mRNA levels of 4 of its predicted target genes. Using RT-qPCR, we detected a 3-fold increase in MYLK expression during GR overexpression (P=0.03). The expression of ACTG2, CNN1 and HLA-DPB1 was not significantly altered between the parental and transfected cells (data not shown).
Discussion
TNBC is a particularly aggressive form of BC, leading to a poor prognosis for patients. Both the absence of hormone receptors and its molecular heterogeneity render TNBC a difficult target for therapeutic intervention. Additionally, a high GR expression was recently linked to therapy failure and worse outcomes in patients with TNBC, as well as other solid tumors. As Chen et al have previously reported, dexamethasone-mediated GR activation induced the expression of genes involved in carcinogenesis and tumor progression (7). These pro-oncogenic processes may also be fostered by GR-responsive non-coding RNAs, including miRNAs. This study therefore aimed at deciphering potential alterations in extracellular and intracellular miRNAs during GR overexpression.
GR biology is fascinatingly complex, involving ligand- dependent receptor activation and isoform-specific transcriptional activity (7,34). In this study, we focused on NR3C1 transcript variant 1, as this is not only the most common isoform in epithelial cells, but also more transcriptionally active than others. Additionally, previous studies have demonstrated notable increases in NR3C1 variant 1 in TNBC (35,36). As expected, the transfection of TNBC cell lines with NR3C1-coding plasmids induced both a strong overexpression of GR mRNA and the increased the expression of several downstream targets (37–39). In concordance with previous reports, we found that TNBC cell lines secrete EVs in an exosome-like size range (15–17). Small RNA from cells with endogenous and artificially increased GR expression was sequenced, and miRNAs compared in differential expression analyses. Although a total of 419 and 442 miRNAs were detected in the MDA-MB-231 and MDA-MB-468 EVs, respectively, there were no statistically significant changes in miRNA profiles for either cell line. Increased metastasis and the progression of GR-overexpressing TNBC may therefore not be mediated by the secretion of soluble mediators to a significant extent. On the other hand, Harris et al recently identified vesicular proteins involved in BC metastasis, indicating a signaling function of vesicle components other than miRNAs (15). Circulating miRNAs are, however, not restricted to EVs, but can also be transported by lipoproteins and protein complexes, such as Argonaute2 (Ago2) (40,41). Given that this study focused exclusively on vesicular miRNAs, no conclusions can be drawn about the impact of GR expression on secreted miRNAs associated with other carrier vehicles.
When assessing the impact of GR overexpression on intracellular expression profiles, we detected a slight, cell line-specific modulation of 7 miRNAs. Even though the impact of GR signaling on TNBC miRNAs has not yet been elucidated, previous studies have reported GR-responsive miRNAs in primary lymphocytes, as well as in liver and spleen cells (42–44). In our data, the non-overlapping profiles of regulated miRNAs in the studied cell lines may be reflective of the inherent transcriptional heterogeneity of TNBC phenotypes (13,45,46). Rainer et al reported similar findings for several lymphoma cell lines that, although featuring GR-reactive miRNAs, only displayed moderate and nonuniform changes in miRNA profiles upon GR activation by dexamethasone (47). In our data, the parental and transfected cells were separated by changes in the global expression profile of miRNAs (Fig. 7). Distances between groups decreased when reducing the number of analyzed miRNAs, indicating that GR slightly shifts the expression patterns of many miRNAs, instead of inducing large changes in the abundance of a few specific transcripts.
Of note, miR-203a-3p, upregulated by GR expression in MDA-MB-436 cells, is controversially discussed in BC literature. Several studies have reported its overexpression in BC, as well as an association with a poor prognosis (48–50). Different data, on the other hand, have suggested that miR-203a-3p serves as a tumor suppressor miRNA, and have stated a decreased expression in BC (51-53). In this study, we found miR-203a-3p to be significantly upregulated in the MDA-MB-436 cells upon the overexpression of GR. As GR is known to be associated with tumor progression, this finding may corroborate the postulation of miR-203a-3p as an oncogenic factor in BC. However, considering the magnitude of expression changes, its biological impact may be of minor relevance. In line with our findings, a previous study reported an upregulation of miR-203 in dexamethasone-treated bone cells, potentially indicating a common miRNA response to GR activation across cell types (54).
The myosin light chain kinase (MYLK, MLCK) has been shown to interact with PI3K-AKT and p38 signaling, increasing cell motility and inhibiting apoptosis in BC cells (55,56). Furthermore, Sundararajan et al pointed out the ability of MYLK to promote invasiveness in several BC cell lines (57). Additionally, using LC-MS/MS-based proteomic profiling, Lawrence et al reported TNBC to feature particularly high levels of MYLK compared to less aggressive BC variants (58). As MYLK is a predicted miR-203a-3p target, we quantified its expression in MDA-MB-436 levels using RT-qPCR. Surprisingly, an increased miR-203a-3p expression in the transfected cells was accompanied by a 3-fold increase in MYLK mRNA levels. This finding was not in concordance with our expectations, as canonical miRNA regulation involves binding of mRNAs and repressing their translation. Consequently, the observed upregulation indicates that MYLK is not directly bound and downregulated by miR-203a-3p.
Taken together, our data suggest that MYLK may be regulated by GR, and can be regarded as a candidate gene involved in the poor survival rates of TNBC patients overexpressing GR. GR's mode of action on MYLK, however, seems not to be mediated by major alterations in cellular miRNAs.
In conclusion, we did not observe any prominent alterations in cellular or vesicular miRNA profiles upon overexpression of GR. The patterns of miRNA expression seem to be influenced by GR to only a small degree, and other mechanisms may therefore be the primary driver for the higher mortality rates of patients suffering from TNBC with GR overexpression.