RNA sequencing identifies specific PIWI-interacting small non-coding RNA expression patterns in breast cancer

Abstract

PIWI-interacting small non-coding RNAs (piRNAs) are genetic and epigenetic regulatory factors in germline cells, where they maintain genome stability, are involved in RNA silencing and regulate gene expression. We found that the piRNA biogenesis and effector pathway are present in human breast cancer (BC) cells and, analyzing smallRNA-Seq data generated from BC cell lines and tumor biopsies, we identified >100 BC piRNAs, including some very abundant and/or differentially expressed in mammary epithelial compared to BC cells, where this was influenced by estrogen or estrogen receptor β, and in cancer respect to normal breast tissues. A search for mRNAs targeted by the BC piRNome revealed that eight piRNAs showing a specific expression pattern in breast tumors target key cancer cell pathways. Evidence of an active piRNA pathway in BC suggests that these small non-coding RNAs do exert transcriptional and post-transcriptional gene regulatory actions also in cancer cells.

Figure 1. Expression of PIWIL proteins and PIWI-interacting RNAs (piRNAs) in BC cell lines

A. Relative abundance of Piwil2 and Piwil4 mRNAs, respect to α-tubulin mRNA, by real-time quantitative rtPCR (top) and of PIWIL2 and PIWIL 4 proteins by Western blot (bottom) in three breast cancer cell lines and in mammary epithelial MCF10A cells. Oligonucleotide sequences and uncropped images can be found in Supplementary Table S1A and Fig. S1. B. Length distribution of unique reads in a representative small RNA sequencing library. C. Visualization of mean-centered and normalized data relative to piRNAs differentially expressed in three breast cancer cell lines respect to mammary epithelial MCF10A cells, with average fold-changes (p-value ≤0.05, Fisher's exact test) shown in the green-red heatmap to the right.

Figure 2. Expression of mRNAs encoding key components of the Piwi/piRNA pathway in MCF-7 cells and modulation of piRNA expression by cell growth and presence of Estrogen Receptor β

A. Relative abundance, respect to α-tubulin mRNA, of transcripts encoding known components of the piRNA biogenesis pathway in wt-MCF7 cells by real-time quantitative rtPCR. Heatmaps showing mean-centered and normalized data relative to piRNAs differentially expressed in exponentially growing (G) vs growth arrested (A) wt MCF7 cells (panel B) in Ct- and Nt-ERβ+ vs wt (ERβ-) MCF-7 cells (panel C) or in exponentially growing vs growth arrested ERβ+ cells (panel D). In all cases, differentially expressed piRNAs identified with a non-parametric Wilcoxon Mann-Whitney test (α <0.05) are shown to the left with hierarchical clustering of replicates by Euclidean distance, together with average fold-changes (Pval ≤0.05, Fisher's exact test) in the green-red heatmaps to the right and the relative principal component analysis (PCA) plots at the bottom.

Figure 3. Identification of 8 piRNAs differently expressed in cancer vs normal breast tissues

A, Heatmaps showing mean-centered and normalized data relative to piRNAs differentially expressed in cancerous (T) vs non-tumoral (N) breast tissue samples from the same patients, identified with a non-parametric Wilcoxon Mann-Whitney test (α <0.05), with average fold-changes (p-value ≤0.05, Fisher's exact test) shown in the green-red heatmap to the right. B, Boxplots summarizing differences in expression of the same piRNAs between tumor and matched non-tumor samples. Samples 40, 45, 53 and 79 correspond to patients TAX577740, TAX577745, TAX577453, and TAX577579, respectively.

Figure 4. Examples of piRNA-RNA sequence complementarity exploited to identify transcripts representing putative targets of the piRNAs found differently expressed in cancer vs normal breast tissue

miRanda was used to identify breast cancer RNAs showing significant sequence complementarity with the eight piRNAs showing significant differences in expression between cancer and matched normal breast tissues. Analysis was performed applying both Watson-Crick ( | ) or GU wobble (: ) base pairing, stringent alignment scores (≥170) and a high binding energy threshold (≤-20.0 kcal/mol). As an example, matches (in red box) for piRNA DQ598677 in the 5′ UTR of TAX1BP mRNA (A), in the 3′ UTR of TNFESF10B mRNA (NM_147187) and in the long non-coding RNA NR_027140 encoded by the same locus (B) and in the coding region of SFRP2 mRNA (C), respectively, are shown. For details, see Supplementary Table S6G.

Figure 5. Functional annotation analysis of biological processes involving mRNAs targeted by 8 piRNAs differentially expressed in cancer vs matched normal breast tissue biopsies

Ingenuity Pathway Analysis (IPA) was used to identify biological processes significantly associated (p-value ≤0.05) with mRNAs targeted by the eight piRNAs found differentially expressed in neoplastic vs normal breast tissue samples from the same patients.