miRNA gene mutations commonly disrupt the proper functioning of miRNA genes

miRNA gene mutations commonly disrupt the proper functioning of miRNA genes | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Article miRNA gene mutations commonly disrupt the proper functioning of miRNA genes Magdalena Machowska, Natalia Szóstak, Adrian Tire, Wladyslaw Wegorek, and 5 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-7029847/v1 This work is licensed under a CC BY 4.0 License Status: Posted Version 1 posted You are reading this latest preprint version Abstract A growing number of mutations are being identified in the noncoding genome, including miRNA genes. However, little is known about the consequences of these mutations and how harmful they are to the functioning of miRNA genes. To evaluate the effects of miRNA gene mutations, we took advantage of a large collection of somatic mutations identified in miRNA genes in > 10K TCGA cancer samples and compared them with the corresponding miRNA-seq data. Using different analytical approaches and highly rigorous statistical criteria, we identified many mutations (n = 87) that affect the level of mature miRNAs (predominantly decreasing), isomiR profiles (precision of DROSHA/DICER1 cleavage), and/or 5p/3p miRNA strand balance. Taken together, the analysis revealed that most miRNA gene mutations, not only those in the seed, may be deleterious for the proper functioning of miRNA genes. We also showed that most miRNA gene mutations destabilize the structure of miRNA precursors and that mutations identified as deleterious are associated with a stronger destabilizing effect. Moreover, although most cancer somatic mutations are randomly occurring neutral variants, some mutations that alter the function of well-known cancer-related miRNA genes, such as MIR21 , MIR142 , or MIR205 , might be functional variants in cancer. Biological sciences/Genetics/Functional genomics Biological sciences/Molecular biology/Non-coding RNAs/miRNAs Biological sciences/Genetics/Gene regulation Biological sciences/Genetics/Mutation Biological sciences/Computational biology and bioinformatics Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Figure 7 Introduction microRNAs (miRNAs) are small (~ 20–23 nt) noncoding regulatory RNAs that are estimated to modulate the expression of almost all human genes 1 . Over three decades of investigations have shown that miRNAs play important roles in the regulation of many cellular and physiological processes, such as cell growth and proliferation, the cell cycle, cell adhesion, apoptosis, cell signaling, and nervous system development 2 , 3 . Numerous miRNAs have been found to be upregulated or downregulated in specific physiological conditions or diseases 4 . The two best-known databases of miRNAs, miRBase 5 and MirGeneDB 6 , cumulatively annotate ~ 2000 miRNAs/miRNA genes in humans, including ~ 600 of high confidence. Although the genetic variation in the noncoding parts of the genome is still largely understudied, an increasing number of genetic variants have been identified in miRNA genes, including common and rare SNPs, sporadic germline mutations, and cancer somatic mutations 7 – 9 . Progress in the identification of noncoding variants in the noncoding genome has been facilitated by the increased popularity of whole-genome sequencing (WGS). The extraction of variants located in miRNA genes from whole-genome datasets and the annotation of these variants may be facilitated by the recently developed miRMut pipeline 10 . However, although there are still few functional analyses of mutations in miRNA genes, those performed have focused almost exclusively on mutations located in the seeds of mature miRNAs, which constitute only a small fraction of miRNA gene/precursor sequences. A few such mutations have been proven to be causative variants of rare Mendelian diseases (germline mutations) 9 , 11 – 14 or to affect the functionality of miRNA genes in cancer (somatic mutations) 15 – 19 . On the other hand, almost nothing is known about the consequences of mutations located in other parts of miRNA genes, especially in their most crucial parts, i.e., the sequence encoding the pre-miRNA hairpin and its immediate flanking sequences. To the best of our knowledge, the only study that has attempted systematic analysis of the genetic variants in miRNA genes was performed 15 years ago (in the pre-NGS era). Dr. Rossi's team, with the use of simple in vitro molecular tests (luciferase assays and northern blots), analyzed a few SNPs identified at that time in different parts of miRNA genes 20 . The results of this study suggested that most of the genetic variants in miRNA genes affect miRNA biogenesis and/or function. Due to the lack of information on the consequences of genetic variants in miRNA genes, distinguishing neutral variants from destructive variants (deleterious to miRNA genes and likely functional) is difficult. This distinction between likely neutral (synonymous) and likely functional or deleterious (missense, nonsense, frameshift splice-site) variants in protein-coding genes greatly facilitates many aspects of research on these genes, including prioritizing variants for further analysis. Canonical miRNAs are generated in a multistep process 21 – 23 . Briefly, miRNA biogenesis begins with a long primary precursor (pri-miRNA) transcribed by RNA polymerase II. The crucial part of each pri-miRNA is a characteristic hairpin structure constituting a secondary precursor (pre-miRNA), which is recognized and excised by the microprocessor complex composed of DGCR8 and the nuclease DROSHA. The pre-miRNA is exported to the cytoplasm, where it is further processed by the nuclease DICER1, which cuts off its terminal loop and generates a miRNA duplex. Once loaded into the miRNA-induced silencing complex (miRISC), the miRNA duplex is unwound, one of its strands is released, and the other becomes the guide strand (mature miRNA). One or both strands of the miRNA duplex may serve as a mature miRNA that, upon complementary interaction with its target sequence, mainly in the 3'UTR of the targeted mRNA, downregulates the expression of the targeted gene by translational repression and/or RNA deadenylation and degradation. Each step of this process may be further regulated by the interaction of miRNA precursors with various RNA-binding proteins. An additional aspect of miRNA biogenesis is that, due to the inaccuracy of DROSHA and/or DICER1 cleavage, more than one miRNA may be generated from some of the miRNA precursors (arms); such miRNAs are termed isomiRs. Although isomiRs may recognize different targets, their functions are mostly unknown, and they are only rarely annotated in miRNA databases (e.g., two isoforms for each arm of miR-142 are annotated in MirGeneDB). Each step of miRNA biogenesis and processing strongly depends on different structural and/or sequence motifs 24 , 25 . The alterations in these motifs via genetic variants/mutations, not only those located in seed regions, may lead to erroneous functioning of miRNA genes, including differences in the efficiency of mature miRNA release and generation of altered miRNAs. The variants may act as riboSNitches, affecting the structure of the precursors 26 , altering the thermodynamic properties of the precursors, and/or directly affecting key functional elements of miRNA precursors, such as DROSHA/DICER1 cleavage sites or the binding motifs of regulatory proteins. Depending on the effect of the variant, it may be considered a loss- or gain-of-function mutation. The further arguments that may support the potential functional/deleterious effect of mutations in miRNA genes are (i) the decreased frequency of SNPs (purifying selection) in these regions 27 , 28 , (ii) the results of testing artificial miRNA precursors/shRNAs, in which the tested sequence changes often resulted in altered processing of the precursors 29 – 32 , and (iii) the results demonstrating the sequence specificity of DROSHA or DICER1 cleavage sites 33 – 35 . To shed more light on the consequences of mutations in miRNA genes (defined as sequences encoding pre-miRNAs with immediately adjacent 25 nt flanking sequences, roughly corresponding/overlapping with miRNA genes designated by the HUGO Gene Nomenclature Committee) on miRNA biogenesis, we took advantage of a large collection (n ~ 7000) of somatic mutations identified and annotated in miRNA genes 7 in cancer samples from The Cancer Genome Atlas (TCGA) project and the corresponding miRNA-seq data annotated with the use of isoMiRmap 36 . IsoMiRmap is a recently developed tool that allows precise annotation of miRNA reads considering length (isomiRs) and sequence (mutations) variation 36 . The comparison of the particular mutations with the corresponding miRNA-seq data allowed us to evaluate in real samples (not artificial functional models) the effect of the substantial number of mutations on miRNA levels, generated isomiRs, and the proportion of generated miRNA strands (5p/3p strand balance). Additionally, with the use of RNA-seq data for the selected mutations, we evaluated the effects of the mutations on the levels of the target mRNAs. The performed analyses revealed that a substantial number of tested mutations severely affect the function of miRNA genes and allow the identification of mutations with the most striking effects on particular aspects of biogenesis. This result is consistent with the idea that miRNA precursors are rather fragile structures in which even subtle changes affect their processing and suggests that most mutations in miRNA genes should be considered potentially functional (deleterious or gain of function) rather than neutral variants. Notably, however, the mutations identified in this study as affecting miRNA genes should not be considered physiologically functional, nor do they necessarily play a role in cancer. Such analyses were beyond the scope of our study. By definition, only a small fraction of mutations occurring in the cancer genome (including mutations in miRNA genes) are expected to play a role in cancer and/or act as cancer-driving mutations. To prove such an effect, a specific set of functional tests focused on and adjusted to the expected impact of a particular mutation (or miRNA) in a specific cancer condition would have to be performed. Methods Data resources and generated metadata The IDs, sequences, and genomic coordinates of miRNAs (5p and 3p strands) and miRNA precursors were obtained from miRBase v.22.1 5 and used as references for mutation annotation and isomiRs classification. The miRNA gene IDs were used according to HUGO Gene Nomenclature. The following cancer type names and abbreviations were used according to TCGA nomenclature: adrenocortical carcinoma (ACC), bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC), cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), lymphoid neoplasm diffuse large B-cell lymphoma (DLBC), esophageal carcinoma (ESCA), head and neck squamous cell carcinoma (HNSC), kidney chromophobe (KICH), kidney renal clear cell carcinoma (KIRC), kidney renal papillary cell carcinoma (KIRP), acute myeloid leukemia (LAML), brain lower grade glioma (LGG), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), mesothelioma (MESO), ovarian serous cystadenocarcinoma (OV), pancreatic adenocarcinoma (PAAD), pheochromocytoma and paraganglioma (PCPG), prostate adenocarcinoma (PRAD), rectum adenocarcinoma (READ), sarcoma (SARC), skin cutaneous melanoma (SKCM), stomach adenocarcinoma (STAD), testicular germ cell tumor (TGCT), thyroid carcinoma (THCA), thymoma (THYM), uterine corpus endometrial carcinoma (UCEC), uterine carcinosarcoma (UCS), and uveal melanoma (UVM). The list of 7110 cancer somatic mutations in miRNA genes, together with sample identifiers and characteristics, was retrieved from Urbanek-Trzeciak et al. (Supplementary Table S2 in 7 ). The mutations were identified in > 10,000 TCGA cancer samples representing 33 cancer types. To analyze the effect of mutations in miRNA genes on the generated miRNAs, we utilized small RNA-seq data of the corresponding TCGA samples mapped with the isoMiRmap tool 36 and retrieved from https://cm.jefferson.edu/isoMiRmap/ 36 . IsoMiRmap allows precise annotation of miRNA length variants (isomiRs) and designated sequence variants, including miRNA reads with non-template nucleotides at the 3' end and miRNA reads with mutations (substitutions and 1-nt indels) annotated in Catalogue of Somatic Mutations in Cancer (COSMIC) v87 (released on 13 November 2018), covering a substantial fraction of the TCGA mutations 7 . For each sample, we downloaded three separate files generated by isoMiRmap, i.e., ‘exclusive’ (containing reads mapped to unique miRNA genes), ‘ambiguous’ (containing reads mapped to multiple genome sites/miRNA genes), and ‘snps’ (containing reads with genetic variants). A detailed description of isoMiRmap, the settings used for mapping, and the format of the output data can be found in 36 . To prevent any uncertainties, biases, or confounding factors, the following mutations were removed from the analysis: (i) mutations identified in redundant miRNA genes (e.g., MIR1-1 (miR-1-1) and MIR1-2 ); (ii) more than one mutation in a particular miRNA gene in one sample; (iii) long indels (more than 4 nt); (iv) mutations in samples with no small RNA-seq data; (v) mutations in miRNA genes with not defined one of the miRNA arms in miRBase; (vi) mutations in samples with more than one small RNA-seq dataset (sequenced multiple times); (vii) mutations in miRNA genes for which no reads were mapped in the mutated samples; (viii) mutations in samples with many ambiguous (mapping in multiple positions) reads in small RNA-seq; (ix) mutations located in mature miRNA or in proximity (± 2 nt) to DROSHA/DICER1 cleavage sites but not annotated in COSMIC v87. On the basis of the analysis of the above data and the above-listed criteria, we created a dataset consisting of 1309 mutations (1218 unique mutations) (Supplementary Table S1). Supplementary Table S1, in addition to data retrieved from 7 (yellowish columns), contains data generated in this study (white column). The new data were generated with the in-house script prepared in R. The list of mutations collected in Supplementary Table S1 served as a base list for selecting mutations for subsequent Experiments (see Fig. 1 and Results). To analyze the effects of mutations on miRNA levels, isomiR distributions, strand balance, and the levels of miRNA targets (mRNAs) (analysis across multiple samples), we used crude small RNA-seq data retrieved from Loher et al. 36 and batch-corrected RNA-seq data retrieved from Hoadley et al. 37 , supplementary file EBPlusPlusAdjustPANCAN_IlluminaHiSeq_RNASeqV2.geneExp.tsv . The small RNA-seq data were batch-corrected by us as previously described 38 , considering all isoMiRmap-annotated isomiRs. We retained all isomiRs to ensure that those potentially unique to mutated samples were not excluded. Batch correction was performed separately for each of the 33 TCGA cancer types (cohorts), considering the following confounding factors: (i) the platform (Illumina Genome Analyzer (GA) or Illumina HiSeq), (ii) tumor purity, (iii) and plate (in which the cDNA library was prepared). We did not correct the data for the type of RNA isolation protocol (direct (total RNA) or MultiMACs (poly-A depleted RNA)) because a single protocol was always used for a specific cancer type. Data about the platform and protocol were collected from the supplementary file in Hoadley et al. 37 and from the GDC Legacy Archive; tumor purity values were estimated using the TCGA tumor purity function from the TCGAbiolinks R package 39 ; and the plate IDs were identified on the basis of the aliquot barcodes. For batch correction, we used two algorithms: ComBat (version 3.80 from the sva R package 40 ) or limma (version 3.6 from the limma R package 41 available as part of the Bioconductor project). To identify potential confounders and select the proper algorithm, principal component analysis (PCA) was performed using the prcomp function of the stats R package ( https://www.rdocumentation.org/packages/stats/versions/3.6.2/topics/prcomp.html ). On the basis of the PCA results, the following statistical tests were performed: the Kruskal‒Wallis test for categorical variables with n > 2, the Wilcoxon rank‒sum test for categorical variables with n = 2, and the Kendall tau rank correlation for ordinal variable purity. The effectiveness of the performed batch corrections was estimated using the FDR p-value, which indicates the association between potential confounders and the principal component (PC), and visualized with the PCA and t-distributed stochastic neighbor embedding (t-SNE) graphs (Rtsne wrapper for the Barnes–Hut t-Distributed Stochastic Neighbor Embedding from the Rtnse R package was used). On the basis of the above, limma was chosen for all tumor cohorts that required correction (CHOL, GBM, and UVM did not show any batch effects and thus were not batch-corrected). miRNA precursor structure modeling The secondary structures of the wild-type and mutant miRNA precursors and the change in Gibbs free energy (revised/optimized dG value) of the structures and individual base pairs or regions were predicted with the use of mfold version 3.6 or the mfold web server 42 with default parameters and processed in VARNA 43 . For modeling, we used pre-miRNA sequences (reconstructed based on miRBase as in 7 , 44 ) extended upstream (5’) and downstream (3’) by 25 nt flanking sequences. The structures with the lowest free energy are presented in this study. Spatial (3D) structures were modeled (on the basis of the predicted secondary structures) using RNAComposer software with the default parameters 45 and visualized via PyMOL 46 . miRNA level analysis Analysis of the influence of miRNA gene mutations on the miRNA level was performed using two different approaches. In the first approach [Experiment_1], within a specific mutant sample, we directly compared the levels (number of reads) of miRNA from the wild-type and the mutant alleles. To evaluate the fold change (FC) (or depletion) of the mutated allele, the numbers of reads at the RNA level were subsequently compared with the corresponding numbers of reads of mutant and wild-type alleles at the DNA level, and the significance of the difference was calculated with Fisher's exact test. This approach was applied only for mutations in mature miRNA sequences that can be distinguished at the RNA level. For the analysis, we considered only the ‘exclusive’ and ’mutated’ type of reads (wild-type and mutant). To select mutations adequate for the test and to avoid any uncertainties, the following mutations were excluded from the analysis: (i) not located in mature miRNA and (ii) located in miRNAs with fewer than 50 total miRNA reads in a mutated sample (not distinguishing reads from wild-type and mutant alleles). In the second approach [Experiment_3], in which we analyzed mutations located in all miRNA gene subregions (not only in mature miRNAs), we compared the level of miRNA in a mutated sample with the levels of miRNA in other samples of the same cancer type (without the mutation). For the analysis, we used batch-corrected RPM values. The level of a specific miRNA gene in a specific sample was calculated as the total RPM of all 5p and 3p isomiRs derived from a particular gene, including ‘exclusive’, ‘non-template’, and ‘mutated’ reads. To increase the normality of the miRNA level distribution, the RPM values were log2 transformed (RPM + 1). To avoid uncertainties, mutations with an RPM value lower than 10, either in a mutated sample or the average RPM of the other samples, were excluded. The level of a miRNA in a sample with a mutation vs. that in a in samples without mutations was compared with a one-sample t-test and z-score. To minimize the chance of false-positive results, only mutations with (i) a t-test with a Bonferroni-corrected p < 0.05, (ii) a z-score p 1, and (iv) identified as outliers defined as a datapoint > 1.5 IQR (interquartile range) below Q1 or above Q3 were considered significant. IsomiR classes Consistently with the convention described previously 47 , all isomiRs were categorized into 9 basic classes denoted as follows: 0|0 (reference (canonical) isomiR as annotated in miRBase); 0|-; 0|+; +|0; -|0; +|-; -|+; +|+; -|-. The signs/values before and after the vertical line correspond to the 5p and 3p miRNA ends (5p|3p) and indicate the direction of the particular end shift, i.e., no change (0), upstream (3p to 5p, (-)), and downstream (+). Additionally, when comparing isomiR profiles between samples (mutant vs. wild-type), reads with non-template 3p-end modifications were considered an additional class of (non-template) isomiRs, denoted (nt) if any modification was present. The extension of the classification may be used to indicate the exact isomiR coordinates, e.g., -2|-1 (indicating an isomiR with the 5p-end shifted by 2 nucleotides upstream and the 3p-end shifted by 1 nucleotide upstream) or by using (n) for any nucleotide change at the 3p-end to emphasize the modifications at the 5p-end. The principles of isomiR classification and denotation are shown in Supplementary Figure S1. IsomiR distribution analysis As in the case of the analysis of the effect of mutations on the level of miRNA, the effect of mutations on isomiR distribution was analyzed using two different approaches. In the first approach, we compared the fraction of isomiRs (isomiR profiles) of mutant and wild-type alleles, which were divided into 9 isomiR classes (as defined above), within one sample. To calculate the number and fraction of reads classified into particular isomiR classes, we used raw read counts (only the ‘exclusive’ and “mutated” type of reads were taken into account). Rare cases of isomiRs/reads for which it could not be determined whether they originated from wild-type or mutant alleles were excluded from the analysis. From the list of mutations (Supplementary Table S1), we excluded mutations that (i) are not located in mature miRNAs and (ii) have < 50 reads derived from either allele. For statistical analysis of changes in isomiR profiles, we used Pearson’s chi-square independence test and Cramer’s V test. For the chi-square test, we added a pseudocount (number of reads + 5) to each isomiR class. Cramer’s V value measures the strength of the association/relationship between two variables. Cramer’s V values 0.4 were interpreted as no association or small (weak), medium, and large (strong) associations, respectively 48 . In the second approach, we compared the isomiRs (isomiR profiles) derived from a mutated gene in the sample with the mutation and in other samples of the same cancer type without the mutation. From the list of mutations (Supplementary Table S1), we excluded mutations that (i) had less than 20 RPM for a particular miRNA strand in the mutated sample and less than 20 RPM for the average in the corresponding wild-type samples. If the criteria were fulfilled for both strands, the isomiR profiles were analyzed for each strand. As described in the subchapter above, isomiRs were divided into 10 classes, and fractions of particular isomiR classes were calculated on the basis of batch-corrected RPM values. For each miRNA of interest, we subsequently calculated the distance (sum of differences of all isomiR fractions; values from 0 to 2 where 0 indicates no difference and 2 indicates the occurrence of completely different isomiRs) between the isomiR profile in a sample with mutation and an average isomiR profile of samples without mutation. The distance calculated for the sample with the mutation was compared (one-sample t-test and z-score) with similarly calculated distances of all samples without the mutation and visualized on a cumulative graph of isomiR profiles. Only mutations with (i) a t-test with a Bonferroni-corrected p < 0.05, (ii) a z-score p 0.2 were considered significant. Strand balance analysis Strand balance analysis was performed by comparing the 5p/3p strand balance (calculated as log2(5p_raw_counts + 1/3p_raw_counts + 1), including ‘exclusive’, ‘non-template’, and ‘mutated’ reads) of a mutated miRNA gene in a mutated sample vs. the strand balance of the gene in other samples of the same cancer type without mutation. To select mutations that were adequate for the test and to avoid any uncertainties, mutations with fewer than 50 total reads for a particular miRNA gene (Supplementary Table S1) were excluded from the analysis. To ensure a reliable estimation of the strand balance in the other (reference) samples, only samples with a total number of reads ≥ 50 were included in the analysis. To identify mutations that significantly affect strand balance, the values were compared with a one-sample t-test and z-score. Only mutations with (i) a t-test with a Bonferroni-corrected p < 0.05, (ii) a z-score p 10% change in strand fractions between mutated sample and average from other WT samples, and (iv) identified as outliers, defined as a datapoint > 1.5 IQR below Q1 or above Q3, were considered significant. Differential expression For differential expression analysis, we used batch-corrected mRNA (RNA-seq) data (see above in section ‘Data resources and generated metadata’). The analysis was performed with the use of the DESeq2-based 49 R pipeline (SARTools) 50 to compare selected samples with mutations affecting the analyzed aspects of miRNA gene biogenesis with other samples of the same cancer type without a mutation. miRNA target prediction miRNA target prediction was performed using TargetScan 51 and mirDB 52 , allowing target prediction for both wild-type and mutant (custom) sequences or isomiRs. For the prediction of mutation-specific targets, we used the most common and/or the most differentiated isomiRs in a mutated sample. We subsequently compared the list of predicted targets with the list of differentially expressed genes between the wild-type and mutant samples. Statistical analysis Statistical analyses were performed using the stats R package or the Real Statistics Resource Pack for MS Excel ( https://real-statistics.com/ ). Results Data processing To evaluate the effects of mutations in miRNA genes on the levels and sequences of generated miRNAs, we compiled data from the list of 7,110 cancer somatic mutations identified in miRNA genes in > 10,000 TCGA cancer samples representing 33 cancer types 7 and from the corresponding miRNA reads (small RNA-seq datasets) mapped using isoMiRmap 36 . IsoMiRmap is a tool that, in addition to mapping wild-type miRNAs, allows reads mapping to mutations identified in TCGA cancer samples, annotated in Catalogue of Somatic Mutations in Cancer (COSMIC) v87 36 (for details, see Methods). To avoid confusion resulting from the ambiguous mapping of miRNA reads, we excluded mutations located in redundant/duplicated miRNA genes [e.g., MIR1-1 (miR-1-1) and MIR1-2 or MIRLET7A1 (let-7a-1), MIRLET7A2 , and MIRLET7A3 ] from the analysis. Additionally, we removed mutations/samples for which there was no small RNA-seq data and mutations in genes for which no miRNA reads were detected in the corresponding samples. Other exclusion criteria are described in the Methods section (“Data resources and generated metadata”) and presented in Fig. 1 . The remaining 1309 mutations listed/characterized in Supplementary Table S1 served as a resource of mutations for subsequent analyses. The effect of mutations on miRNAs expressed from the mutated genes was analyzed using two general approaches, i.e., comparing miRNAs expressed from a particular miRNA gene from the mutant and wild-type alleles in the mutated sample (APPROACH_1) and comparing miRNAs expressed from a mutated gene in the sample with the mutation with other samples of the same cancer type without the mutation (APPROACH_2). The general concept of the study, as well as the exclusion criteria and the number of mutations selected for subsequent Experiments, are shown in Fig. 1 and Supplementary Figure S2. Analysis of the effect of miRNA gene mutations by direct comparison of miRNAs generated from mutated and wild-type alleles in individual samples with mutations In this approach, we took advantage of the fact that some mutations may be directly observed at the RNA level; therefore, direct observation of the miRNAs derived from the mutant and wild-type alleles in one sample is possible without biases resulting from the biological and technical differences between samples or their processing/analysis. The limitation of such an approach, however, is that it allows the analysis only of mutations located in sequences encoding mature miRNAs (expressed at the miRNA level). In this approach, to ensure proper distinction between wild-type and mutant reads, only mutations annotated in COSMIC v87 were taken into account. In the first experiment (Experiment_1), to investigate the effects of mutations on miRNA levels, we compared the allelic fractions of the mutations (proportions of reads with mutations) at the genomic (DNA) and transcript (mature miRNA) levels. To ensure proper determination of the mutant fraction, only mutations in miRNAs with at least 50 miRNA reads (wild-type + mutant) were used for analysis (n = 53). Surprisingly, as shown in Fig. 2 A, there was very little relationship between the fractions of mutant alleles at the DNA and RNA levels, and most mutations (n = 32, 60%) significantly (adjusted p < 0.05; Fisher’s exact test) deviated from the trendline (x = y), representing an equal proportion of mutant reads at the DNA and miRNA levels, expected for mutations with no effect on the miRNA level. Only one mutation significantly increased the miRNA level (adj. p > 0.05; absolute FC ≥ 2), whereas the vast majority of these mutations decreased the miRNA level of mutated alleles (1 vs. 31; p = 0.00003; Fisher’s exact test), in most cases almost to 0. Moreover, the relative decreases (FCs) had much larger amplitudes than the relative increases in miRNA levels (Fig. 2 A inset). As shown in Fig. 2 B, the mutations affecting miRNA levels are roughly equally distributed along mature miRNA sequences and do not cluster in any specific position of the precursor (details regarding each of the mutations are listed in Supplementary Table S2). It is noteworthy that several mutations (n = 5) that occurred in more than one sample, especially in samples of the same cancer type, generally had consistent effects on miRNA levels (Fig. 2 A, inset), further confirming the validity of the analysis. Among these mutations, n.22G > T in MIR379 was identified in two LUAD samples and, in both cases, acutely decreased the level of mutant miRNA, whereas n.75C > T in MIR342 (in 2 SKCM samples), n.59T > C in MIR142 (in 2 DLBC samples), and n.51G > A in MIR320A (in LUAD, LUSC, and HNSC samples) had no significant effect on the miRNA level. Only n.35C > T in MIR205 had a discordant effect on the miRNA level in the SKCM sample (decrease) vs. the CESC and LUSC samples (no significant change). These discrepancies may, however, result from differences in cancer/tissue type, differences in mutation allelic frequency, and general variation in cancer samples. Recently, it was shown that the nucleotide composition of a 3-base pair-long motif of the miRNA precursor duplex encompassing the DICER1 cleavage site, called GYM, affects the efficiency of DICER1 processing and thus may affect miRNA levels 53 . By comparing the GYM scores of the 3 mutations located in the GYM region with those of their WT counterparts, we found that the effect of one of the mutations, i.e., n.35C > T in MIR523 , which decreases the miRNA level, may be well explained by the GYM effect, i.e., the GYM score decreased from 26 for the WT precursor to 10 for the mutant precursor (the GYM score ranges from 0-100, with higher values indicating more efficient processing; Fig. 2 C). Notably, however, miRNA levels can be influenced by many other factors, including the structure of the miRNA precursor (discussed below). In Experiment_2, to investigate the impact of mutations on the precision of DROSHA/DICER1 cleavage, i.e., generated isomiRs, we compared the isomiR profiles of particular miRNAs generated from wild-type and mutant alleles. To ensure the reliable determination of isomiR distributions, for the analysis, we selected only mutations/samples with at least 50 miRNA reads for both the wild-type and mutant alleles (n = 16; Supplementary Table S3). To compare the distribution of the corresponding wild-type and mutant miRNAs, we classified all the miRNA reads into 9 isomiR classes categorized on the basis of the position of their ends (5p|3p) upstream (+) or downstream (-), against the corresponding ends annotated in miRBase (0|0) (as proposed previously 47 and graphically illustrated in Supplementary Figure S1). Ten of the 16 (63%) tested mutations induced significant (chi-square adj. p 0.2, indicating at least a moderate relationship of the mutations with the isomiR profile (Supplementary Table S3). Among the identified mutations, 3 indicate a strong effect on the isomiR profile (Cramer’s V > 0.4). The mutations affecting isomiR profiles are roughly equally distributed along mature miRNA sequences (Fig. 3 A). To determine whether Cramer’s V value reliably distinguishes changes induced by mutations from random variation in isomiR distributions, we also compared isomiR profiles of the same mutants and wild-type alleles of the same miRNAs in different samples. As shown in Fig. 3 B, the Cramer’s V values are much lower for pairs of the same mutants or wild-type alleles than for pairs of corresponding mutant and wild-type alleles from the same samples, even though the former come from different samples or even different cancer types. Thus, a Cramer's V threshold of 0.2 accurately distinguishes significant changes from random variability. Most of the identified mutations induce changes of > 2-fold in at least one isomiR class and, in most cases, substantially (> 20%) change the level of the main isomiRs (Fig. 3 C, D, E; Supplementary Figure S3). Notably, 6 of the mutations induced changes predominantly at the 5p-end, resulting in a shift in the miRNA seed sequences, which have a direct impact on target recognition. Among the mutations with the strongest effect on the isomiR distribution was n.75C > T in MIR342 , located at the 15th nucleotide of miR-342-3p, which shortened the miRNA at the 3p-end and thus severely reduced (> 4-fold) the fraction of the canonical (0|0) and 0|+ isomiRs in favor of the 0|- isomiR (Fig. 3 C). This mutation was identified in two SKCM samples. The very similar isomiR profile in two independent samples with this mutation (Cramer’s V > 0.2) proves the validity of the results and the specificity of the changes induced by the mutation. An example of a mutation affecting the 5p-end of the generated miRNAs is n.55A > G in MIR142 , located at the 3rd nucleotide of the miR-142-3p seed (Cramer’s V = 0.638). This mutation elongates the 5p-end of the generated isomiRs by one nucleotide, severely reducing the level of the +|- and +|0 isomiRs (predominant in the wild-type allele) and increasing the level of the 0|- and 0|0 isomiRs (Fig. 3 D). A similar shift (in favor of isomiRs that are longer at the 5p end) is induced by 3 other mutations located at nearby positions in the miR-142-3p seed (n.58G > C and n.59T > C (2x; in 2 samples), although the strength of the effect varies depending on the mutation and is not formally significant (Cramer’s V ~ 0.18) in all cases (Fig. 3 E). We subsequently combined all the isomiRs of the MIR142 mutant samples on the basis of particular changes at the 5p-end (5p-isomiRs; ignoring changes at the 3p-end). As shown in Fig. 3 F, the predominant isomiR class expressed from the wild-type allele (accounting for approximately 80% of all 4 samples) consists of isomiRs + 1|n shortened by 1 nt at the 5p-end. The + 1|n isomiRs are expressed at much lower levels in the mutant alleles, ranging from 9 to 61%, depending on the mutation. In contrast, the mutant isomiRs consisted of a greater fraction (38–89%) of isomiRs with a canonical 5p-end (0|n), accounting for only ~ 20% of the wild-type isomiRs (Fig. 3 F). This shift of the 5p-miRNA end (+ 1|n to 0|n) results in an additional change in the seed sequence that is independent of the point changes directly introduced by the individual mutations. The remaining mutations that significantly affect isomiR distribution are shown in Supplementary Figure S3. Analysis of the effect of miRNA gene mutations by comparing miRNAs in samples with the mutation versus corresponding samples without the mutation To extend the analysis to mutations located in other parts of miRNA genes (not only in mature miRNAs that can be observed at the RNA level), we used another approach, in which we analyzed the effect of mutations by comparing miRNAs expressed from a mutated gene in a sample with the mutation with other samples of the same cancer type that lack the mutation. However, as the power of such an approach is strongly limited by the dilution of the mutation effect by the presence of the normal allele and the contamination of cancer samples with normal (non-cancerous) cells (in most cases, a mutant allele accounts for < < 50%; Supplementary Table S1), as well as very high genetic and transcriptional variation in cancer samples, we sacrificed sensitivity to focus on identifying single mutations with the most profound effects. For this purpose, we used the above standard multiple criteria for classifying results as significant (see Methods), minimizing the potential of false-positive results but likely missing the effects of many mutations (false-negatives). In Experiment_3, to identify mutations that affect the level of miRNAs, we compared batch-corrected levels of miRNAs (log2-transformed RPM values of the 5p and 3p arms) in samples with mutations to the average levels of the miRNAs in samples without mutations. As shown in Fig. 4 A, of the 682 mutations (648 unique mutations) that fulfilled the criteria of the analysis, 21 significantly affected the miRNA level, including 15 mutations that increased the level and 6 mutations that decreased the level (Fig. 4 , Supplementary Table S4). The observed excess of mutations increasing the level of miRNA (compared to Experiment_1) is likely attributable to the reduced statistical power to detect downregulation relative to upregulation. As shown in Fig. 4 B, mutations affecting miRNA levels are located in all subregions of the miRNA precursors and are roughly equally distributed along the sequence. Among the examples of mutations affecting miRNA levels are mutation n.5-7delAGC in the 5p-flank in MIR122 , which causes a decrease in the miRNA level (log2(FC) A, which is located in the 3p-flank in MIR518E , which leads to a dramatic increase in the miRNA level (log2(FC) > 10) (Fig. 4 D). Both mutations affect the predicted 2D structure and stability of the miRNA precursors by relaxing the duplex structure near the mutation site. Moreover, both mutations alter the expression of numerous genes, including the predicted targets of the relevant miRNAs (Fig. 4 E, F; Supplementary Table S5). In Experiment_4, to investigate the effect of mutations on the precision of DROSHA/DICER1 cleavage, for each mutation, we compared the isomiR profiles of a given miRNA in the sample with the mutation and corresponding (same cancer type) samples without mutations. To ensure the reliable determination of isomiR profiles, only samples with at least 20 RPM were considered in the analysis. Among the 420 mutations (397 unique mutations) that fulfilled the criteria, 32 mutations significantly affected the isomiR profile in at least 1 strand (Fig. 5 A, Supplementary Table S6). Mutations affecting the isomiR profile are located in all subregions of the miRNA precursors and are roughly equally distributed along the sequence (Fig. 5 B). Among the mutations affecting the isomiR profile is n.85G > A, which is located in the 3p-flank of MIR142 (Fig. 5 C). The mutation leads to extension by one nucleotide of the 5p-end of the generated miR-142-3p isomiRs, severely reducing the level of the +|0 and +|- isomiRs (predominant in the wild-type samples) and increasing the level of the 0|0 and 0|- isomiRs. Interestingly, other mutations located in different parts of MIR142 induced a similar shift in the isomiR profile (as analyzed in Experiment_2, compare the isomiR profile graph in Figure ). Another mutation that alters the isomiR profile is n.21-22delTA in MIR10B , which is located in the 5p flank of the gene (Fig. 5 D). The mutation causes the extension of the 5p-end of miR-10b-5p by one or two nucleotides, leading to a reduction in the fraction of 0|n isomiRs in favor of + 1|n and + 2|n isomiRs. Another example is n.-6C > G in MIR21 , which is located in the 5p flank of the gene (Fig. 5 E). The mutation induces shortening by 1 nt of the 3p miRNA at its 3p-end, severely reducing the fraction of 0|+ isomiRs predominant in WT samples in favor of the canonical (0|0) isomiR. As shown in the volcano plots (Supplementary Figure S4), some genes are differentially expressed in the mutated samples. However, the relationship between the differentially expressed genes and the observed isomiR shifts cannot be directly determined. The isomiR profiles of other mutations that significantly affect isomiR profiles are presented in Supplementary Figure S5. In Experiment_5, we analyzed the effect of mutations on miRNA strand balance, i.e., the ratio of miRNAs derived from 5p and 3p miRNA arms of precursor, compared between mutant and wild-type samples as log2(5p + 1/3p + 1). To ensure reliable determination of the strand balance, only samples with ≥ 50 reads of a given miRNA were considered. As shown in Fig. 6 A, from the 515 mutations (490 unique mutations) fulfilling the criteria of the analysis (for details, see Methods), 4 mutations significantly affected the miRNA strand balance, including one mutation located in each of the following regions: 5p-flank, 3p-flank, loop, and miRNA duplex (Fig. 6 A, B; Supplementary Table S7). The most striking example is the duplication of two Gs (n.97_98dupGG) in the 3p-flank in MIR205 identified in THYM. In all the THYM samples, miR-205-5p was the predominant strand, accounting for more than 99.7% of all the reads (similar to other TCGA cancers, as well as other tissues reported in different databases). In contrast, in the mutant sample, the fraction of miR-205-5p decreased substantially (down to 65%), resulting in a dramatic increase (> 100x, up to 35%) of the miR-205-3p fraction (Fig. 6 C). It is noteworthy that the observed change of the 5p/3p balance is detectable even though the fraction of the mutant allele (at the DNA level) accounts for only 19%, which suggests that the actual effect of the mutant allele is much stronger (this comment also applies to other mutations mentioned in this section). A closer look at the miRNAs generated from the 5p and 3p arms revealed that, in addition to changing the 5p/3p balance, the mutation also shifted (by 1 to 3 nucleotides toward the base of the precursor hairpin) the predominant DROSHA/DICER1 cleavage sites, generating abnormal miRNAs (Fig. 6 C). To determine whether this dramatic strand/isoform imbalance (the mutant-specific miRNA) affects gene expression, we performed differential expression analysis and identified 31 and 5 genes whose expression was downregulated and upregulated, respectively, in the sample with the MIR205 mutation. Among the downregulated genes, four targets (predicted by TargetScan and/or mirDB) of mutant-specific miR-205-5p (-2|-1) and 3p (+ 2|+3) were identified (Supplementary Table S5; Fig. 6 E). An example of a mutation that disturbs the balance in the opposite direction (in favor of the 5p strand) is n.84C > T in MIR365B , which is located in the 3p-strand of the miRNA duplex (Fig. 6 D). This mutation increases the level of the 5p-strand more than fivefold, from 3–19%. As the mutation destabilizes the miRNA duplex at the 5p-end of miR-365b-5p (replaces G:C with the G:U pair and increases the duplex free energy by 2 kcal/mol), its effect of favoring the 5p over the 3p strand is consistent with the principle of greater preference/stability (more efficient loading into the miRISC complex) of the strand with a more relaxed (less stable) 5p-end 54 (Fig. 6 D). As shown in Fig. 6 F, many genes were downregulated in the mutated sample, and 7 of these genes were predicted targets of miR-365b-5p, the fraction of which was increased in the mutant sample. The other two mutations that affect the strand balance are shown in Supplementary Figure S6. The impact of miRNA gene mutations on the structure of miRNA precursors Although mutations can affect the functioning of miRNA genes through various mechanisms, the most direct/common impact seems to be modifying the structure of miRNA precursors. Therefore, to preliminarily estimate the effects of mutations on the structure of miRNA precursors, we modeled and compared the wild-type and mutant precursor structures of all the tested miRNA genes (n = 1309). Most mutations increase the change in the Gibbs free energy (dG) value between mutant and WT precursors, i.e., decrease the stability of the precursor structures (Supplementary Table S8). The average dG value of the mutant structures was 0.77 kcal/mol greater than that of the corresponding wild-type structures (p = 1.5E-35; paired t-test), and there was an excess of mutations increasing the dG value (destabilizing the precursor structure) over mutations decreasing the dG value (Fig. 7 A). Further comparisons revealed that the increase in dG (ddG) was significantly greater for mutations affecting the functioning of miRNA genes (functional mutations identified in Experiments_1–5) than for the remaining mutations (t-test, p = 1.73E-08). As shown in Fig. 7 B, the functional mutations cluster at the top of the ddG ranking. The ddG for functional mutations is independent of the functional effect of the mutation; it applies both to mutations impacting miRNA levels (t-test, p = 6.92E-09) and mutations affecting isomiR distributions (t-test, p = 0.0012) (Fig. 7 B; due to the low number of mutations, we did not test mutations affecting 5p/3p-strand balance). The relationship between the results of structural and functional analyses further confirms the reliability of the identified functional mutations (results of Experiments_1–5), as both analyses are completely independent. Examples of mutations affecting the miRNA precursor structures are shown in Fig. 7 C. Discussion Thanks to the genetic code, predicting the consequences of mutations and distinguishing deleterious from neutral mutations in coding sequences is relatively straightforward (although not trivial). On the other hand, owing to modern sequencing approaches, such as WGS, thousands of new genetic variants, primarily located in noncoding sequences, are being identified in each sequenced sample. Unfortunately, due to the heterogeneity of noncoding elements and a lack of appropriate "codes" (rules), predicting the consequences of mutations or even estimating the likelihood of deleteriousness for mutations in the best-defined noncoding elements, such as miRNA genes, is nearly impossible. Consequently, the vast majority of sequence variants identified in the noncoding genome are assumed to be neutral or variants of unknown significance (VUS) and are ignored. To overcome these limitations and ultimately develop "codes" to predict the consequences of variants in noncoding elements, more data need to be collected on the impacts of mutations on the functionality of these elements. Among the most well-recognized and highly conserved noncoding elements are miRNA genes, particularly their ~ 100 nt long regions encoding hairpin-structured miRNA precursors (pre-miRNAs). Although analysis of genetic variation in the noncoding genome has long been neglected, numerous variants in miRNA genes have been extracted from whole-genome datasets. Moreover, few miRNA gene mutations have been implicated in Mendelian diseases, including mutations in MIR96 in non-syndromic hearing loss 11 , mutations in MIR184 in different hereditary eye diseases 12 , mutations in MIR204 in retinal dystrophy 13 , and mutations in MIR140 in skeletal dysplasia 14 . There are also a few miRNA genes whose mutations have been suggested to be potential cancer drivers, including the MIR15A/MIR16-1 cluster in CLL, MIR142 in different blood cancers, MIR122 in liver carcinoma, and MIR21 in various cancers 7 , 15 , 55 – 58 . Nevertheless, a larger-scale analysis of the molecular consequences of miRNA gene variants has never been performed, and almost all previous studies focused on mutations identified in the seed sequences 11 , 13 , 14 , 17 , 18 , 59 – 62 . Only a few studies have investigated the impact of variants in other regions of miRNAs 19 , 20 , 25 , 63 – 67 (summarized in 9 ). Moreover, almost all of these studies have investigated only the impact on the miRNA level. To date, the only and largest systematic analysis of sequence variants in miRNA genes, performed more than 15 years ago, has examined the effects of 24 common SNPs located in different parts of miRNA precursors. Using several cellular functional assays, the study revealed that most of the variants disturb the proper function of miRNA genes, including the efficiency of miRNA generation and precision of miRNA processing 20 . Therefore, taking advantage of the large collection of mutations and corresponding transcriptome (miRNA-seq) datasets generated by the TCGA project, we performed the most comprehensive analysis to date of the impact of mutations on miRNA gene functioning, exceeding the number of variants studied in previous projects by orders of magnitude. In contrast to most previous studies, the approach allowed us, for the first time, to analyze the effects of mutations in their real genetic context, under natural conditions in real cells/tissues, without artificially generated models using engineered genes, overexpression, or generic cell lines. Furthermore, in contrast to most previous studies, our analysis included mutations located in all parts of the miRNA gene, not only in the seed region. An additional advantage of our study is that the vast majority of somatic mutations are randomly occurring variants. This ensures that the analysis is not biased toward more likely neutral variants, such as common SNPs, or deleterious variants, such as mutations identified in Mendelian diseases. Moreover, in addition to the analysis of the impact of mutations on the miRNA level, which has been mostly analyzed in other studies, we also analyzed the effects of mutations on the miRNA strand balance and isomiR profiles. Altogether, we analyzed the effects of 703 mutations selected in multiple steps from ~ 7000 mutations identified in miRNA genes in TCGA samples (Fig. 1 , Supplementary Figure S2). The mutations were selected on the basis of very strict criteria, ensuring that their effects may be reliably evaluated (e.g., by excluding mutations in redundant miRNA genes such as MIR1-1 and MIR1-2 ). Nevertheless, this is the largest and most complex analysis of miRNA gene mutations performed to date. First, by comparing the levels of mutant and wild-type alleles at the genomic (DNA) and transcriptomic (miRNA) levels, we showed that most miRNA gene mutations (at least those located in the miRNA duplex) affect the miRNA level (32/53) and/or precision of DROSHA/DICER1 cleavages (isomiR profile 10/16), i.e., they are deleterious for the proper functioning of miRNA genes at the molecular level. Additional mutations, which are located in the whole pre-miRNA and affect miRNA gene functionalities, including miRNA levels (n = 21), isomiR profiles (n = 32), and miRNA strand balance (n = 4), were detected by comparing the level and sequence of miRNAs in the mutated vs. non-mutated samples. However, due to the dilution of mutants with the wild-type allele, contamination of cancer samples with normal tissues, and cancer heterogeneity, the sensitivity of the latter approach was low. Therefore, the mutations identified here likely reflect only a fraction of the functionally relevant events and should be interpreted as representative mutations rather than an exhaustive set. In none of the analyses did we observe that mutations affecting the miRNA genes clustered in any subregion of the gene, which could indicate the particular importance of this region for the tested functionalities. However, interesting observations may include the mutation n.84C > T in MIR365B , which changes the stability of one side of miRNA duplexes and thus drastically reverses the miRNA strand balance. This observation is in agreement with the notion that the miRNA strand with the less stable 5p-end is more favorably loaded into the RISC, becoming a more stable mature miRNA 54 . Finally, we showed that most miRNA gene mutations decrease the stability of miRNA precursors. We also showed that mutations that we detected as having a functional effect cause, on average, a greater decrease in miRNA precursor stability. This association suggests that at least some mutation effects are expressed via changes in the miRNA precursor structure and confirms the importance of the structure and its stability for the proper functioning of miRNA precursors. Although most tested mutations, including those affecting the molecular functionality of miRNA genes, are randomly occurring neutral variants, single ones may still play a role in driving cancer development. This applies, in particular, to mutations in cancer-related miRNA genes, such as those annotated in the Cancer miRNA Census 55 , including 28 miRNA genes annotated in CMC, in which mutations affecting miRNA gene functioning were detected. One possible example of such a gene is MIR142 , which is recurrently mutated in various hematological malignancies, including AML, CLL, DLBCL, FL, and other types of B-cell lymphomas, and was found as the most frequently mutated miRNA gene in any cancer 7 , 68 . MIR142 is highly expressed in lymphoid blood cells, and it has been shown that miR-142 (particularly miR-142-3p) plays an important role in hematopoiesis, regulating the development and function of different hematologic lineages as well as in hematological malignancies 69 . Functional studies have focused mainly on a few mutations in the miR-142-3p seed. However, these studies, conducted in various cellular and animal models, demonstrated that the mutations, through ineffective regulation of miR-142-3p targets, including ASH1L , increase the level of HOXA9/A10, resulting in aberrant hematopoietic differentiation 17 . This promotes myeloid and suppresses lymphoid lineages, ultimately leading to leukemic transformation and AML. Mutations synergize with mutations in IDH2 59 . A gain-of-function effect for one of the mutations was also suggested 70 . The functional effects of the seed mutations were observed despite their small impact on the level of miR-142-3p, suggesting that the mutations act predominantly by affecting the seed sequence and thus target recognition. This finding is consistent with the lack of effects of n.59T > C (two samples) and n.58G > C on the miR-142-3p level observed in our study (Experiment_1; Fig. 2 A). A recent study published during the preparation of this manuscript also demonstrated the molecular effects of mutations located outside the 3p seed in different parts of MIR142 overexpressed in HEK293 cells 19 . To date, however, MIR142 mutations have been studied only in artificially generated cellular or mouse models. Here, for the first time, we demonstrated the molecular effects of MIR142 mutations in relevant cancer samples, predominantly LAML and DLBC. We found that two mutations, n.55A > G and n.16C > A, severely affected the levels of miR-142-3p and miR-142-5p (Experiment_1; Supplementary Table S2). Additionally, three MIR142 mutations significantly affected the distribution of miR-142-3p isomiRs (n.55A > G and n.59T > C in Experiment_2 and n.85G > A in Experiment_4). Surprisingly, despite being located in different parts of the precursor, all these mutations induced similar effects on the isomiR profile, leading to a decrease in + 1|n isomiRs (predominant in wild-type) in favor of 0|n isomiRs. A similar (although less profound) shift in the isomiR profile was also triggered by two other mutations in MIR142 (n.58G > C and n.59T > C). Thus, all of these mutations affect the miR-142-3p seed by shifting the miRNA; additionally, 3 of the mutations located in the 3p seed also directly alter its sequence. Another well-known cancer-related miRNA gene is MIR205 , which is important for cancer development. It has been assigned both an oncogenic role and a tumor suppressor role, depending on the tissue/cancer type 71 , 72 . Analysis of TCGA datasets revealed enrichment of mutations in MIR205 7 in different solid cancers, particularly in melanoma; however, their impact on miRNA functioning has never been studied. Moreover, a reduced level of miR-205-5p is associated with shorter survival in melanoma patients 73 . In this study, we found that n.35C > T, located in the 5p-arm of MIR205 (Experiment_1), decreases the level of miR-205 in SKCM (melanoma). Two other mutations, n.30C > T and n.94G > C, identified in BLCA and localized in 5p- and 3p-flanks, respectively, affect the isomiR profiles of miR-205-5p (Experiment_4; Supplementary Figure S5). Finally, n.97-98dupGG, located in the 3p-flank (Experiment_5, Fig. 6 C), affects the structure of miRNA precursors, resulting in a shift in DROSHA and DICER1 cleavage sites and reversing the 5p/3p strand balance, reducing the fraction of miR-205-5p (dominant mature strand) in favor of miR-205-3p (passenger strand). Two more examples of mutations located in cancer-related genes are n.21-22delTA in MIR10B 74 and n.-6C > G in MIR21 75 . Both mutations led to significant changes in isomiR profiles (Experiment_4; Fig. 5 ). Mutation in MIR10B increases the + 1|n and + 2|n fractions of miR-10b-5p, changing its canonical seed. In contrast, the mutation in MIR21 affects miR-21-3p, decreasing the fraction of 0|+ in favor of 0|0 and non-templated isomiRs. In summary, to the best of our knowledge, we performed the broadest systematic analysis of the effects of mutations in miRNA genes, in which we identified 87 mutations that significantly affect different miRNA gene functionalities, including miRNA levels, isomiR profiles, and strand balance. For the first time, we studied the effects of mutations in natural genomic contexts, not in artificially generated models. Our results improve the understanding of the impact of genetic variants on miRNA biogenesis and may help develop tools for predicting the significance of genetic variants in miRNA genes. Generally, our results (the first approach) indicate that most miRNA gene mutations, not only those located in seeds, affect the proper functioning of miRNA genes and should therefore be considered likely to reveal deleterious variants in genetic analyses. Declarations Authors’ contribution Magdalena Machowska : Conceptualization, Methodology, Validation, Formal analysis, Resources, Data Curation, Visualization, Writing - Original Draft, Writing - Review & Editing. Natalia Szostak : Methodology, Software, Validation, Formal analysis, Resources, Data Curation, Writing - Review & Editing. Adrian Tire : Formal analysis, Visualization, Writing - Review & Editing. Wladyslaw Wegorek : Formal analysis, Visualization, Writing - Review & Editing. Malwina Suszynska : Formal analysis, Writing - Review & Editing. Arkadiusz Kajdasz : Software, Formal analysis, Writing - Review & Editing. Paulina Galka-Marciniak : Formal analysis, Funding acquisition, Writing - Review & Editing. Anna Philips : Resources, Writing - Review & Editing. Piotr Kozlowski : Conceptualization, Methodology, Validation, Formal analysis, Resources, Data Curation, Visualization, Supervision, Project administration, Funding acquisition, Writing - Original Draft, Writing - Review & Editing. Acknowledgements The study was supported by grants from the Polish National Science Centre (2020/39/B/NZ5/01970 and 2016/22/A/NZ2/00184 awarded to Piotr Kozlowski, and 2020/39/D/NZ2/03106 awarded to Paulina Galka-Marciniak). Competing interests The authors have no conflicts of interest to declare. References Gebert LFR, MacRae IJ (2019) Regulation of microRNA function in animals. Nat Rev Mol Cell Biol 20:21–37 O’Brien J, Hayder H, Zayed Y, Peng C (2018) Overview of MicroRNA Biogenesis, Mechanisms of Actions, and Circulation. Front Endocrinol (Lausanne) 9 Saliminejad K, Khorshid K, Soleymani Fard HR, S., Ghaffari SH (2019) An overview of microRNAs: Biology, functions, therapeutics, and analysis methods. J Cell Physiol 234:5451–5465 Ali Syeda Z, Langden SSS, Munkhzul C, Lee M, Song SJ (2020) Regulatory Mechanism of MicroRNA Expression in Cancer. Int J Mol Sci 21:1723 Kozomara A, Birgaoanu M, Griffiths-Jones S (2019) miRBase: from microRNA sequences to function. Nucleic Acids Res 47:D155–D162 Fromm B et al (2022) MirGeneDB 2.1: toward a complete sampling of all major animal phyla. Nucleic Acids Res 50:D204–D210 Urbanek-Trzeciak MO et al (2020) Pan-cancer analysis of somatic mutations in miRNA genes. EBioMedicine 61:103051 Han M, Zheng Y (2013) Comprehensive Analysis of Single Nucleotide Polymorphisms in Human MicroRNAs. PLoS ONE 8:e78028 Machowska M, Galka-Marciniak P, Kozlowski P (2022) Consequences of genetic variants in miRNA genes. Comput Struct Biotechnol J 20:6443–6457 Urbanek-Trzeciak MO, Kozlowski P, Galka-Marciniak P, miRMut (2022) Annotation of mutations in miRNA genes from human whole-exome or whole-genome sequencing. STAR Protoc 3:101023 Mencía Á et al (2009) Mutations in the seed region of human miR-96 are responsible for nonsyndromic progressive hearing loss. Nat Genet 41:609–613 Hughes AE et al (2011) Mutation Altering the miR-184 Seed Region Causes Familial Keratoconus with Cataract. Am J Hum Genet 89:628–633 Conte I et al (2015) MiR-204 is responsible for inherited retinal dystrophy associated with ocular coloboma. Proceedings of the National Academy of Sciences 112 Grigelioniene G et al (2019) Gain-of-function mutation of microRNA-140 in human skeletal dysplasia. Nat Med 25:583–590 Calin GA et al (2002) Frequent deletions and down-regulation of micro- RNA genes miR15 and miR16 at 13q14 in chronic lymphocytic leukemia. Proceedings of the National Academy of Sciences 99, 15524–15529 Cimmino A et al (2005) miR-15 and miR-16 induce apoptosis by targeting BCL2. Proceedings of the National Academy of Sciences 102, 13944–13949 Trissal MC et al (2018) MIR142 Loss-of-Function Mutations Derepress ASH1L to Increase HOXA Gene Expression and Promote Leukemogenesis. Cancer Res 78:3510–3521 Kwanhian W et al (2012) Micro RNA -142 is mutated in about 20% of diffuse large B‐cell lymphoma. Cancer Med 1:141–155 Bayraktar R et al (2025) The mutational landscape and functional effects of noncoding ultraconserved elements in human cancers. Sci Adv 11 Sun G et al (2009) SNPs in human miRNA genes affect biogenesis and function. RNA 15:1640–1651 Leitão AL, Enguita FJ (2022) A Structural View of miRNA Biogenesis and Function. Noncoding RNA 8:10 Treiber T, Treiber N, Meister G (2019) Regulation of microRNA biogenesis and its crosstalk with other cellular pathways. Nat Rev Mol Cell Biol 20:5–20 Shang R, Lee S, Senavirathne G, Lai EC (2023) microRNAs in action: biogenesis, function and regulation. Nat Rev Genet 24:816–833 Baek SC et al (2024) Structural atlas of human primary microRNAs generated by SHAPE-MaP. Mol Cell 84:1158–1172e6 Roden C et al (2017) Novel determinants of mammalian primary microRNA processing revealed by systematic evaluation of hairpin-containing transcripts and human genetic variation. Genome Res 27:374–384 Halvorsen M, Martin JS, Broadaway S, Laederach A (2010) Disease-associated mutations that alter the RNA structural ensemble. PLoS Genet 6:e1001074 Omariba G et al (2020) Genome-Wide Analysis of MicroRNA-related Single Nucleotide Polymorphisms (SNPs) in Mouse Genome. Sci Rep 10:5789 Gong J et al (2012) Genome-wide identification of SNPs in microRNA genes and the SNP effects on microRNA target binding and biogenesis. Hum Mutat 33:254–263 Fang W, Bartel DP (2015) The Menu of Features that Define Primary MicroRNAs and Enable De Novo Design of MicroRNA Genes. Mol Cell 60:131–145 Rice GM, Shivashankar V, Ma EJ, Baryza JL, Nutiu R (2020) Functional Atlas of Primary miRNA Maturation by the Microprocessor. Mol Cell 80:892–902e4 Auyeung VC, Ulitsky I, McGeary SE, Bartel DP (2013) Beyond secondary structure: primary-sequence determinants license pri-miRNA hairpins for processing. Cell 152:844–858 Galka-Marciniak P et al (2016) siRNA release from pri-miRNA scaffolds is controlled by the sequence and structure of RNA. Biochim Biophys Acta 1859:639–649 Starega-Roslan J, Witkos T, Galka-Marciniak P, Krzyzosiak W (2015) Sequence Features of Drosha and Dicer Cleavage Sites Affect the Complexity of IsomiRs. Int J Mol Sci 16:8110–8127 Nguyen TD, Trinh TA, Bao S, Nguyen TA (2022) Secondary structure RNA elements control the cleavage activity of DICER. Nat Commun 13:2138 Kim K et al (2021) A quantitative map of human primary microRNA processing sites. Mol Cell 81:3422–3439e11 Loher P et al (2021) IsoMiRmap: fast, deterministic and exhaustive mining of isomiRs from short RNA-seq datasets. Bioinformatics 37:1828–1838 Hoadley KA et al (2018) Cell-of-Origin Patterns Dominate the Molecular Classification of 10,000 Tumors from 33 Types of Cancer. Cell 173:291–304e6 Ibing S et al (2021) On the impact of batch effect correction in TCGA isomiR expression data. NAR Cancer 3 Colaprico A et al (2016) TCGAbiolinks: an R/Bioconductor package for integrative analysis of TCGA data. Nucleic Acids Res 44:e71–e71 Leek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD (2012) The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics 28:882–883 Ritchie ME et al (2015) limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res 43:e47–e47 Zuker M (2003) Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res 31:3406–3415 Darty K, Denise A, Ponty YVARNA (2009) Interactive drawing and editing of the RNA secondary structure. Bioinformatics 25:1974–1975 Galka-Marciniak P et al (2019) Somatic Mutations in miRNA Genes in Lung Cancer—Potential Functional Consequences of Non-Coding Sequence Variants. Cancers (Basel) 11:793 Popenda M et al (2012) Automated 3D structure composition for large RNAs. Nucleic Acids Res 40:e112–e112 The PyMOL Molecular Graphics System, Version 3.0 Schrödinger, LLC. Loher P, Londin ER, Rigoutsos I (2014) IsomiR expression profiles in human lymphoblastoid cell lines exhibit population and gender dependencies. Oncotarget 5:8790–8802 Rea LM, Parker RA (1992) Designing and Conducting Survey Research. Jossey-Bass, San Francisco Love MI, Huber W, Anders S (2014) Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 15:550 Varet H, Brillet-Guéguen L, Coppée J-Y, Dillies M-A, SARTools (2016) A DESeq2- and EdgeR-Based R Pipeline for Comprehensive Differential Analysis of RNA-Seq Data. PLoS ONE 11:e0157022 Agarwal V, Bell GW, Nam J-W (2015) & Bartel, D. P. Predicting effective microRNA target sites in mammalian mRNAs. Elife 4 Chen Y, Wang X (2020) miRDB: an online database for prediction of functional microRNA targets. Nucleic Acids Res 48:D127–D131 Lee Y-Y, Kim H, Kim V (2023) N. Sequence determinant of small RNA production by DICER. Nature 615:323–330 Medley JC, Panzade G, Zinovyeva A (2021) Y. microRNA strand selection: Unwinding the rules. WIREs RNA 12 Suszynska M et al (2024) Cancer miRNA Census – a list of cancer-related miRNA genes. Nucleic Acids Res 52:1628–1644 Dietlein F et al (2022) Genome-wide analysis of somatic noncoding mutation patterns in cancer. Sci (1979) 376 Hornshøj H et al (2018) Pan-cancer screen for mutations in non-coding elements with conservation and cancer specificity reveals correlations with expression and survival. NPJ Genom Med 3:1 Rheinbay E et al (2020) Analyses of non-coding somatic drivers in 2,658 cancer whole genomes. Nature 578:102–111 Marshall A et al (2020) Mir142 loss unlocks IDH2R140-dependent leukemogenesis through antagonistic regulation of HOX genes. Sci Rep 10:19390 Soldà G et al (2012) A novel mutation within the MIR96 gene causes non-syndromic inherited hearing loss in an Italian family by altering pre-miRNA processing. Hum Mol Genet 21:577–585 Duan R, Pak C, Jin P (2007) Single nucleotide polymorphism associated with mature miR-125a alters the processing of pri-miRNA. Hum Mol Genet 16:1124–1131 Ding W et al (2018) A polymorphism rs3746444 within the pre-miR‐499 alters the maturation of miR‐499‐5p and its antiapoptotic function. J Cell Mol Med 22:5418–5428 Ohanian M, Humphreys DT, Anderson E, Preiss T, Fatkin D (2013) A heterozygous variant in the human cardiac miR-133 gene, MIR133A2, alters miRNA duplex processing and strand abundance. BMC Genet 14:18 Zhao D et al (2022) Exosomal miR-1304-3p promotes breast cancer progression in African Americans by activating cancer-associated adipocytes. Nat Commun 13:7734 Lechner J et al (2013) Mutational Analysis of MIR184 in Sporadic Keratoconus and Myopia. Invest Opthalmology Visual Sci 54:5266 Calin GA et al (2005) A MicroRNA Signature Associated with Prognosis and Progression in Chronic Lymphocytic Leukemia. N Engl J Med 353:1793–1801 Iwai N, Naraba H (2005) Polymorphisms in human pre-miRNAs. Biochem Biophys Res Commun 331:1439–1444 Galka-Marciniak P et al (2022) Mutations in the miR-142 gene are not common in myeloproliferative neoplasms. Sci Rep 12:10924 Huang W, Paul D, Calin GA, Bayraktar R (2023) miR-142: A Master Regulator in Hematological Malignancies and Therapeutic Opportunities. Cells 13:84 Kawano S et al (2023) A gain-of‐function mutation in microRNA 142 is sufficient to cause the development of T‐cell leukemia in mice. Cancer Sci 114:2821–2834 Chauhan N, Dhasmana A, Jaggi M, Chauhan SC, Yallapu MM (2020) miR-205: A Potential Biomedicine for Cancer Therapy. Cells 9:1957 Ferrari E, Gandellini P (2020) Unveiling the ups and downs of miR-205 in physiology and cancer: transcriptional and post-transcriptional mechanisms. Cell Death Dis 11:980 Sánchez-Sendra B et al (2018) Downregulation of intratumoral expression of miR-205, miR-200c and miR-125b in primary human cutaneous melanomas predicts shorter survival. Sci Rep 8:17076 Sheedy P, Medarova Z (2018) The fundamental role of miR-10b in metastatic cancer. Am J Cancer Res 8:1674–1688 Bautista-Sánchez D et al (2020) The Promising Role of miR-21 as a Cancer Biomarker and Its Importance in RNA-Based Therapeutics. Mol Ther Nucleic Acids 20:409–420 Additional Declarations There is NO Competing Interest. Magdalena Machowska and Natalia Szostak are first joint authors. Supplementary Files MachowskaMSupplementaryFigures02072025.pdf Supplementary Figures S1-S6 MachowskaMSupplementaryTableS102072025.xlsx Supplementary Table S1 MachowskaMSupplementaryTableS202072025.xlsx Supplementary Table S2 MachowskaMSupplementaryTableS302072025.xlsx Supplementary Table S3 MachowskaMSupplementaryTableS402072025.xlsx Supplementary Table S4 MachowskaMSupplementaryTableS502072025.xlsx Supplementary Table S5 MachowskaMSupplementaryTableS602072025.xlsx Supplementary Table S6 MachowskaMSupplementaryTableS702072025.xlsx Supplementary Table S7 MachowskaMSupplementaryTableS802072025.xlsx Supplementary Table S8 Cite Share Download PDF Status: Posted Version 1 posted You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-7029847","acceptedTermsAndConditions":true,"allowDirectSubmit":true,"archivedVersions":[],"articleType":"Article","associatedPublications":[],"authors":[{"id":485051096,"identity":"1ebdbf1c-c8a8-46ce-95f3-2a80d4354620","order_by":0,"name":"Magdalena Machowska","email":"","orcid":"https://orcid.org/0000-0003-0118-6241","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland; 2. Faculty of Biotechnology, Univesrity of Wroclaw, Wroclaw, Poland","correspondingAuthor":false,"prefix":"","firstName":"Magdalena","middleName":"","lastName":"Machowska","suffix":""},{"id":485051097,"identity":"05458edc-629c-4feb-9af2-b9d70e2ae95e","order_by":1,"name":"Natalia Szóstak","email":"","orcid":"https://orcid.org/0000-0001-9933-0496","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Natalia","middleName":"","lastName":"Szóstak","suffix":""},{"id":485051098,"identity":"83051d02-8aee-47ce-9e69-88e01e9969ad","order_by":2,"name":"Adrian Tire","email":"","orcid":"","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Adrian","middleName":"","lastName":"Tire","suffix":""},{"id":485051099,"identity":"283efda6-3e54-4454-b5d9-5c53f6409883","order_by":3,"name":"Wladyslaw Wegorek","email":"","orcid":"","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Wladyslaw","middleName":"","lastName":"Wegorek","suffix":""},{"id":485051100,"identity":"a72cf4d2-6a8c-4ea4-a6e8-3a38430a67c0","order_by":4,"name":"Malwina Suszynska","email":"","orcid":"","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Malwina","middleName":"","lastName":"Suszynska","suffix":""},{"id":485051101,"identity":"aff5cbf1-d63d-4e07-a52d-db895ca62dd7","order_by":5,"name":"Arkadiusz Kajdasz","email":"","orcid":"","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Arkadiusz","middleName":"","lastName":"Kajdasz","suffix":""},{"id":485051102,"identity":"406a118b-98c9-4f25-8223-d6cbca777e31","order_by":6,"name":"Paulina Galka-Marciniak","email":"","orcid":"https://orcid.org/0000-0003-0729-8746","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Paulina","middleName":"","lastName":"Galka-Marciniak","suffix":""},{"id":485051103,"identity":"8acdb428-412f-4727-b7cd-9c4ee4ac5832","order_by":7,"name":"Anna Philips","email":"","orcid":"https://orcid.org/0000-0001-5020-8836","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":false,"prefix":"","firstName":"Anna","middleName":"","lastName":"Philips","suffix":""},{"id":485051095,"identity":"a35d782c-e6a1-42f5-8183-99599f2e56ef","order_by":8,"name":"Piotr Kozlowski","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAA70lEQVRIiWNgGAWjYBACCQYGxgNAWo5NgoENSbwArxYGkBZjNC0GhLUkNhCtRbK9/cGBj2130vukm589YNxhY2/OfvYAcwEeLdI8BxIOzmx7ltsmc8zcgPFMGrNlT14C8ww8WuQkEg4c5t12OLdNIodNgrHtMJvBgRwDZh68WhIbQFrS2aBaeAzOv8GvRVoimQGkJQGmRcLgBgFbJHuOMRyc+e+ZIdAvZhKJZ9IMDG68MTiMzy8Sx9sfPvhw5o68/OzmZxIfgSFmcD7H8HFBBW4tUHAAQiU2QOjDBDXAtTBCtTAToWUUjIJRMApGDgAA8TtPzLxB9moAAAAASUVORK5CYII=","orcid":"https://orcid.org/0000-0003-3770-7715","institution":"Institute of Bioorganic Chemistry, Polish Academy of Sciences, Poznan, Poland","correspondingAuthor":true,"prefix":"","firstName":"Piotr","middleName":"","lastName":"Kozlowski","suffix":""}],"badges":[],"createdAt":"2025-07-02 13:31:01","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-7029847/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-7029847/v1","draftVersion":[],"editorialEvents":[],"editorialNote":"","failedWorkflow":false,"files":[{"id":87177090,"identity":"b473f685-323e-4f8c-be7f-aa00f3fe187d","added_by":"auto","created_at":"2025-07-21 08:54:42","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":349235,"visible":true,"origin":"","legend":"\u003cp\u003eWorkflow scheme and mutation-selection criteria for the analysis.\u003c/p\u003e","description":"","filename":"1.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/ada88880705e61301735a3b7.png"},{"id":87177469,"identity":"75ced4a5-c09a-4253-b742-dee43c643048","added_by":"auto","created_at":"2025-07-21 09:02:48","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":335759,"visible":true,"origin":"","legend":"\u003cp\u003eEffects of miRNA gene mutations on miRNA levels (Experiment_1). A – Comparison of the proportion of mutant reads at the DNA and miRNA levels; inset – volcano plot showing log2(FC) values of the fraction of mutated miRNAs compared with the fraction of mutated allele at the DNA level (x-axis) and -log10(adj.p value) of Fisher's exact test (y-axis). Each dot represents one mutation; red, green, and gray dots indicate mutations increasing, decreasing, and not significantly changing the miRNA level, respectively. The color circles indicate mutations of the same type occurring in different samples. B – Distribution of the analyzed mutations along the consensus miRNA precursor structure. C – The structure of the \u003cem\u003eMIR523\u003c/em\u003e precursor with the GYM motif altered by the n.35C\u0026gt;T mutation.\u003c/p\u003e","description":"","filename":"2.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/d789cb0c091a26a893798e71.png"},{"id":87177112,"identity":"7412a9a8-42a4-43aa-a01c-12c9574849a0","added_by":"auto","created_at":"2025-07-21 08:54:46","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":335882,"visible":true,"origin":"","legend":"\u003cp\u003eEffects of miRNA gene mutations on isomiR profiles (Experiment_2). A – Distribution of the analyzed mutations along the consensus miRNA precursor structure. B – Box and whisker plots showing the distribution of Cramer's V values for differences in isomiR profiles between mutant (MUT) and wild-type (WT) alleles from the same samples, WT alleles of the same miRNAs from different samples, and MUT alleles of the same type from different samples. C – Comparison of WT and n.75C\u0026gt;T miR-342-3p isomiR profiles in 2 SKCM samples. Above: 2D structure of the \u003cem\u003eMIR342\u003c/em\u003e precursor with analyzed mutations indicated with arrowheads; below: miR-342-3p isomiRs representing particular isomiR classes (indicated on the right); on the right: fraction of particular isomiR classes of the WT (black bars indicating an average of two samples; error bars indicating upper and lower values) and the n.75C\u0026gt;T mutant alleles (red and green bars from SKCM_1 and SKCM_2 samples, respectively); on the right: log2(FC) of a particular isomiR class fraction of the mutant alleles against the WT alleles. D – Comparison of WT and n.55A\u0026gt;G miR-142-3p isomiR profiles in an LAML sample. Above: 2D structure of the \u003cem\u003eMIR142\u003c/em\u003eprecursor with the indicated positions of the 4 mutations; below: miR-142-3p isomiRs representing particular isomiR classes (indicated on the right); on the right: fraction of particular isomiR classes of the WT (black bars) and the n.55A\u0026gt;G mutant alleles (red bars); and on the right: log2(FC) of particular isomiR class fractions of the mutant allele against the WT allele. E – Log2(FC) values of isomiR classes of 4 mutations (indicated in Panel D; against WT) located in the miR-142-3p seed. F – Fractions of isomiR classes grouped only on the basis of the 5p-miRNA end (5p-isomiRs). The first bar represents the average value (from 4 samples) for the WT allele, and the following bars represent four \u003cem\u003eMIR142\u003c/em\u003e mutant alleles. The error bar on the WT line indicates extreme values of the isomiR classes observed for the WT allele in four samples with mutations. On the right: sequences of isomiRs included in the 0|n and +|n classes.\u003c/p\u003e","description":"","filename":"3.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/58ad262344f970579b2534b8.png"},{"id":87177089,"identity":"c822a302-9302-446c-a24d-eb98c60e0fed","added_by":"auto","created_at":"2025-07-21 08:54:42","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":384409,"visible":true,"origin":"","legend":"\u003cp\u003eEffects of miRNA gene mutations on miRNA levels (Experiment_3). A – Volcano plot showing log2(FC) values of the miRNA level in samples with mutations in comparison to the average miRNA level in samples (of the same cancer type) without the mutation (x-axis) against -log10(p) of the z-score. Each dot represents a mutation; red, green, and gray (open and closed) dots represent mutations that cause increases, decreases, and no significant change, respectively; the closed gray dots indicate mutations that do not meet the criterion of being outliers in the miRNA level distribution (see Methods). B – Distribution of the analyzed mutations along the consensus miRNA precursor structure. C – Comparison of the \u003cem\u003eMIR122\u003c/em\u003e level in the sample with the n.5-7delAGC mutation with other LIHC samples; from the top: (i) 2D structures of the \u003cem\u003eMIR122\u003c/em\u003e WT and mutant precursors with the indicated position of the n.5_7delAGC mutation; (ii) the juxtaposed graph comparing the dG values of the corresponding base pairs/structural motifs of the WT and mutant precursors; and (iii) the graph showing the distribution of the miR-122 (5p+3p) level in the LIHC samples (pink and blue bars indicate the sample with the mutation and wild-type samples, respectively), inset – box and whisker plot showing log2-transformed levels of the mutant and wild-type samples. D – Comparison of \u003cem\u003eMIR518E\u003c/em\u003e levels in the sample with the n.85G\u0026gt;A mutation vs. other LUAD samples (the panel scheme as in C). E – Volcano plot illustrating differential expression analysis of a single sample with n.5_7delAGC in \u003cem\u003eMIR122\u003c/em\u003evs. corresponding wild-type LIHC samples; red, green, and gray dots indicate upregulated, downregulated, and not significantly changed genes, respectively. The predicted direct targets of miR-122-5p are indicated on the graph in black circles. F – Differential expression analysis of a sample with n.85G\u0026gt;A in \u003cem\u003eMIR518E\u003c/em\u003e vs. corresponding wild-type LUAD samples; the predicted direct targets of miR-518e-5p and 3p are indicated on the graph in black circles (others as in E).\u003c/p\u003e","description":"","filename":"4.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/395a8b12d1d754bbd5121fde.png"},{"id":87177125,"identity":"41b6c374-15f6-4dc3-8731-d187baa095a1","added_by":"auto","created_at":"2025-07-21 08:54:48","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":344911,"visible":true,"origin":"","legend":"\u003cp\u003eEffects of miRNA gene mutations on isomiR profiles (Experiment_4). A – Scatter plot showing the distances between the isomiR profiles of the mutants and their corresponding wild-type samples (x-axis) against the -log10(p) of the z-scores of these distances. Each dot represents a mutation; red and gray dots represent mutations that significantly affect the isomiR profiles and not causing significant changes, respectively. B – Distribution of the analyzed mutations along the consensus miRNA precursor structure. C – Comparison of WT and n.85G\u0026gt;A miR-142-3p isomiR profiles in DLBC samples. From the left: 2D structure of the \u003cem\u003eMIR142\u003c/em\u003e precursor with the indicated position of the mutation and juxtaposed miR-142-3p isomiRs assigned to particular isomiR classes (indicated on the right); a fraction of particular isomiR classes of the WT (average in blue, individual samples in grayscale) and the n.85G\u0026gt;A mutant (pink); log2(FC) of particular isomiR class fractions of the mutant sample against the average of WT samples; a distance of the mutant and WT isomiR profiles from the average WT profile; and fractions of isomiR classes grouped only based on the 5p-miRNA end (5p-isomiRs). D – Comparison of WT and n.21-22delTA miR-10b-5p isomiR profiles in COAD samples (panel scheme as in C). E – Comparison of WT and n.-6C\u0026gt;G miR-21-3p isomiR profiles in DLBC samples (panel scheme as in C except for the lack of a chart showing fractions of 5p-end isomiRs, as the mutation affects only the 3p-end).\u003c/p\u003e","description":"","filename":"5.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/919210a592d9fe0f077ddb82.png"},{"id":87177133,"identity":"4729a976-16d4-406c-80c3-772893926c69","added_by":"auto","created_at":"2025-07-21 08:54:51","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":359427,"visible":true,"origin":"","legend":"\u003cp\u003eEffects of miRNA gene mutations on the balance of miRNA strands (Experiment_5). A – A scatterplot showing the relationship between the effect of mutation on the change in the 5p/3p strand fraction between mutant and wild-type samples (x-axis) and the -log10(p) of the z-score for strand balance (log2(5p+1/3p+1)) (y-axis) in a sample with a mutation compared with other samples of the same cancer type without mutation. Each dot represents a mutation; red and gray (open and closed) dots indicate mutations affecting 5p/3p strand balance and not inducing significant changes, respectively; the closed gray dots indicate mutations not meeting the criterion of being outliers in the strand balance distribution (see Methods). B – Distribution of the analyzed mutations along the consensus miRNA precursor structure. C – Comparison of 5p/3p strand balance in the sample with the n.97-98dupGG mutation in \u003cem\u003eMIR205\u003c/em\u003e vs. other THYM samples without mutation; from the top: (i) 2D structures of the \u003cem\u003eMIR205\u003c/em\u003e WT and mutant precursors with the indicated position of the n.97_98dupGG mutation; (ii) the juxtaposed graph comparing dG values of corresponding wt and mutant base pairs/structural motifs; (iii) the graph showing the proportion of miR-205-5p and miR-205-3p strands in the mutant sample (pink) and in other (wild-type) samples of THYM (blue); inset – box and whisker plot showing log2-transformed 5p/3p strand balance in the mutant sample vs. distribution of the balance in other THYM samples. D – Comparison of 5p/3p strand balance in the sample with the n.84C\u0026gt;T mutation in \u003cem\u003eMIR365B\u003c/em\u003evs. other SKCM samples (the panel scheme as in C). E – Volcano plots illustrating differential expression analysis of a single sample with n.97_98dupGG in \u003cem\u003eMIR205\u003c/em\u003e vs. corresponding wild-type THYM samples (red, green, and gray dots indicate upregulated, downregulated, and not significantly changed genes, respectively). The direct targets of hsa-mir-205-5p -2|n and hsa-mir-205-3p +2|n isomiRs (mostly profound in a mutant sample) are indicated in black circles. F – Differential expression analysis of a sample with n.84C\u0026gt;T in \u003cem\u003eMIR365B\u003c/em\u003e vs. corresponding wild-type SKCM samples (scheme of the panel as in E). The direct targets of the hsa-mir-365b-5p 0|n isomiR are indicated by black circles.\u003c/p\u003e","description":"","filename":"6.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/29dd24458c2a92ed6d1aabee.png"},{"id":87177113,"identity":"0291ea31-baeb-488e-8f20-43a785071ffe","added_by":"auto","created_at":"2025-07-21 08:54:46","extension":"png","order_by":7,"title":"Figure 7","display":"","copyAsset":false,"role":"figure","size":442616,"visible":true,"origin":"","legend":"\u003cp\u003eThe impact of miRNA gene mutations on the stability and structure of miRNA precursors. A – A chart showing the effect of mutations on the stability of miRNA precursor structures, i.e., the change in dG between the mutated and the corresponding wild-type miRNA precursors, expressed as ddG. The graph shows all mutations preselected for the analysis (n=1309; see Figure 1) sorted from the highest (destabilizing mutations) to lowest (stabilizing mutations) ddG values. Light and dark blue lollipops represent mutations significantly affecting the functioning of miRNA genes in one or more Experiments, respectively. In the inset, flowcharts separately highlight only mutations affecting the level of miRNA (orange lollipops, Experiments_1 and 3) and isomiR profiles (green lollipops, Experiments_2 and 4), respectively. B – Box plots comparing the effect of mutations on miRNA precursor stability, i.e., the distribution of ddG values ​​for mutations for which no effect on miRNA gene function was detected (gray box), mutations for which an effect was detected in at least one Experiment (blue box), mutations affecting miRNA levels (orange box), and mutations affecting the isomiR profiles (green box). C – Effect of representative mutations (indicated in A) affecting the functioning of miRNA genes on the secondary (above) and spatial (below) structures of miRNA precursors. Wild-type and mutant structures are shown in dark gray and green, respectively; mutation positions are indicated in red.\u003c/p\u003e","description":"","filename":"7.png","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/a005dd12248533cd6dc25ed5.png"},{"id":88281784,"identity":"56f7d47b-a5e8-46d3-b261-c880f10367a0","added_by":"auto","created_at":"2025-08-04 20:35:14","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":2943172,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/10d8e2ec-0cbc-454b-ba14-9f1fb23b1f81.pdf"},{"id":87177145,"identity":"a6119d17-7c6b-4a45-a840-498790baa689","added_by":"auto","created_at":"2025-07-21 08:54:54","extension":"pdf","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":1136999,"visible":true,"origin":"","legend":"Supplementary Figures S1-S6","description":"","filename":"MachowskaMSupplementaryFigures02072025.pdf","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/cb6867385fc517b50f8cc2f1.pdf"},{"id":87177080,"identity":"0e507c0a-adf9-4c8c-a31f-27305d98eb1a","added_by":"auto","created_at":"2025-07-21 08:54:39","extension":"xlsx","order_by":2,"title":"","display":"","copyAsset":false,"role":"supplement","size":131262,"visible":true,"origin":"","legend":"\u003cp\u003eSupplementary Table S1\u003c/p\u003e","description":"","filename":"MachowskaMSupplementaryTableS102072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/fb7d438fefe102e807a03e9d.xlsx"},{"id":87177118,"identity":"d6491b81-05bb-4cf3-ad7f-4d8541ea6117","added_by":"auto","created_at":"2025-07-21 08:54:47","extension":"xlsx","order_by":3,"title":"","display":"","copyAsset":false,"role":"supplement","size":22972,"visible":true,"origin":"","legend":"Supplementary Table S2","description":"","filename":"MachowskaMSupplementaryTableS202072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/758e93d8ba90e044af1630f5.xlsx"},{"id":87177115,"identity":"d47bbb10-9df9-43e6-b080-773e2ad46800","added_by":"auto","created_at":"2025-07-21 08:54:46","extension":"xlsx","order_by":4,"title":"","display":"","copyAsset":false,"role":"supplement","size":12889,"visible":true,"origin":"","legend":"Supplementary Table S3","description":"","filename":"MachowskaMSupplementaryTableS302072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/6867cc3ea7c07f6dfa9d5754.xlsx"},{"id":87177127,"identity":"3aead16b-26b7-4996-a4f2-6b20fc1fbb3f","added_by":"auto","created_at":"2025-07-21 08:54:48","extension":"xlsx","order_by":5,"title":"","display":"","copyAsset":false,"role":"supplement","size":186116,"visible":true,"origin":"","legend":"Supplementary Table S4","description":"","filename":"MachowskaMSupplementaryTableS402072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/1b57cad1e5a233f46744dca8.xlsx"},{"id":87177091,"identity":"8e0eda1f-4d58-4f37-af2f-139332779415","added_by":"auto","created_at":"2025-07-21 08:54:43","extension":"xlsx","order_by":6,"title":"","display":"","copyAsset":false,"role":"supplement","size":5129825,"visible":true,"origin":"","legend":"Supplementary Table S5","description":"","filename":"MachowskaMSupplementaryTableS502072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/8ecb88fa990f9bca45f2597e.xlsx"},{"id":87177146,"identity":"ec698f37-7f83-4cf2-9604-1f2456a67914","added_by":"auto","created_at":"2025-07-21 08:54:54","extension":"xlsx","order_by":7,"title":"","display":"","copyAsset":false,"role":"supplement","size":146665,"visible":true,"origin":"","legend":"Supplementary Table S6","description":"","filename":"MachowskaMSupplementaryTableS602072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/5579b98f4b370b03dcc6bcc2.xlsx"},{"id":87177111,"identity":"8b75ad29-b8ea-4b17-a995-21139257a4c8","added_by":"auto","created_at":"2025-07-21 08:54:46","extension":"xlsx","order_by":8,"title":"","display":"","copyAsset":false,"role":"supplement","size":185121,"visible":true,"origin":"","legend":"Supplementary Table S7","description":"","filename":"MachowskaMSupplementaryTableS702072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/92954547a191abfddca770b5.xlsx"},{"id":87177083,"identity":"ee3549f5-cf7b-4885-af76-2c62364d9b43","added_by":"auto","created_at":"2025-07-21 08:54:40","extension":"xlsx","order_by":9,"title":"","display":"","copyAsset":false,"role":"supplement","size":153684,"visible":true,"origin":"","legend":"Supplementary Table S8","description":"","filename":"MachowskaMSupplementaryTableS802072025.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-7029847/v1/095aa1678b4a9fd68bbb605f.xlsx"}],"financialInterests":"\u003cp\u003eThere is \u003cstrong\u003eNO\u003c/strong\u003e Competing Interest.\u003c/p\u003e\n\u003cp\u003eMagdalena Machowska and Natalia Szostak are first joint authors.\u003c/p\u003e","formattedTitle":"miRNA gene mutations commonly disrupt the proper functioning of miRNA genes","fulltext":[{"header":"Introduction","content":"\u003cp\u003emicroRNAs (miRNAs) are small (~ 20–23 nt) noncoding regulatory RNAs that are estimated to modulate the expression of almost all human genes\u003csup\u003e\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e\u003c/sup\u003e. Over three decades of investigations have shown that miRNAs play important roles in the regulation of many cellular and physiological processes, such as cell growth and proliferation, the cell cycle, cell adhesion, apoptosis, cell signaling, and nervous system development\u003csup\u003e\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e,\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e\u003c/sup\u003e. Numerous miRNAs have been found to be upregulated or downregulated in specific physiological conditions or diseases\u003csup\u003e\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e\u003c/sup\u003e. The two best-known databases of miRNAs, miRBase\u003csup\u003e\u003cspan citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e\u003c/sup\u003e and MirGeneDB\u003csup\u003e\u003cspan citationid=\"CR6\" class=\"CitationRef\"\u003e6\u003c/span\u003e\u003c/sup\u003e, cumulatively annotate ~ 2000 miRNAs/miRNA genes in humans, including ~ 600 of high confidence.\u003c/p\u003e \u003cp\u003eAlthough the genetic variation in the noncoding parts of the genome is still largely understudied, an increasing number of genetic variants have been identified in miRNA genes, including common and rare SNPs, sporadic germline mutations, and cancer somatic mutations\u003csup\u003e\u003cspan additionalcitationids=\"CR8\" citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e–\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e\u003c/sup\u003e. Progress in the identification of noncoding variants in the noncoding genome has been facilitated by the increased popularity of whole-genome sequencing (WGS). The extraction of variants located in miRNA genes from whole-genome datasets and the annotation of these variants may be facilitated by the recently developed miRMut pipeline\u003csup\u003e\u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eHowever, although there are still few functional analyses of mutations in miRNA genes, those performed have focused almost exclusively on mutations located in the seeds of mature miRNAs, which constitute only a small fraction of miRNA gene/precursor sequences. A few such mutations have been proven to be causative variants of rare Mendelian diseases (germline mutations)\u003csup\u003e\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e,\u003cspan additionalcitationids=\"CR12 CR13\" citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e–\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e\u003c/sup\u003e or to affect the functionality of miRNA genes in cancer (somatic mutations)\u003csup\u003e\u003cspan additionalcitationids=\"CR16 CR17 CR18\" citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e–\u003cspan citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e\u003c/sup\u003e. On the other hand, almost nothing is known about the consequences of mutations located in other parts of miRNA genes, especially in their most crucial parts, i.e., the sequence encoding the pre-miRNA hairpin and its immediate flanking sequences. To the best of our knowledge, the only study that has attempted systematic analysis of the genetic variants in miRNA genes was performed 15 years ago (in the pre-NGS era). Dr. Rossi's team, with the use of simple \u003cem\u003ein vitro\u003c/em\u003e molecular tests (luciferase assays and northern blots), analyzed a few SNPs identified at that time in different parts of miRNA genes\u003csup\u003e\u003cspan citationid=\"CR20\" class=\"CitationRef\"\u003e20\u003c/span\u003e\u003c/sup\u003e. The results of this study suggested that most of the genetic variants in miRNA genes affect miRNA biogenesis and/or function. Due to the lack of information on the consequences of genetic variants in miRNA genes, distinguishing neutral variants from destructive variants (deleterious to miRNA genes and likely functional) is difficult. This distinction between likely neutral (synonymous) and likely functional or deleterious (missense, nonsense, frameshift splice-site) variants in protein-coding genes greatly facilitates many aspects of research on these genes, including prioritizing variants for further analysis.\u003c/p\u003e \u003cp\u003eCanonical miRNAs are generated in a multistep process\u003csup\u003e\u003cspan additionalcitationids=\"CR22\" citationid=\"CR21\" class=\"CitationRef\"\u003e21\u003c/span\u003e–\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e\u003c/sup\u003e. Briefly, miRNA biogenesis begins with a long primary precursor (pri-miRNA) transcribed by RNA polymerase II. The crucial part of each pri-miRNA is a characteristic hairpin structure constituting a secondary precursor (pre-miRNA), which is recognized and excised by the microprocessor complex composed of DGCR8 and the nuclease DROSHA. The pre-miRNA is exported to the cytoplasm, where it is further processed by the nuclease DICER1, which cuts off its terminal loop and generates a miRNA duplex. Once loaded into the miRNA-induced silencing complex (miRISC), the miRNA duplex is unwound, one of its strands is released, and the other becomes the guide strand (mature miRNA). One or both strands of the miRNA duplex may serve as a mature miRNA that, upon complementary interaction with its target sequence, mainly in the 3'UTR of the targeted mRNA, downregulates the expression of the targeted gene by translational repression and/or RNA deadenylation and degradation. Each step of this process may be further regulated by the interaction of miRNA precursors with various RNA-binding proteins. An additional aspect of miRNA biogenesis is that, due to the inaccuracy of DROSHA and/or DICER1 cleavage, more than one miRNA may be generated from some of the miRNA precursors (arms); such miRNAs are termed isomiRs. Although isomiRs may recognize different targets, their functions are mostly unknown, and they are only rarely annotated in miRNA databases (e.g., two isoforms for each arm of miR-142 are annotated in MirGeneDB).\u003c/p\u003e \u003cp\u003eEach step of miRNA biogenesis and processing strongly depends on different structural and/or sequence motifs\u003csup\u003e\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e,\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e\u003c/sup\u003e. The alterations in these motifs via genetic variants/mutations, not only those located in seed regions, may lead to erroneous functioning of miRNA genes, including differences in the efficiency of mature miRNA release and generation of altered miRNAs. The variants may act as riboSNitches, affecting the structure of the precursors\u003csup\u003e\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e\u003c/sup\u003e, altering the thermodynamic properties of the precursors, and/or directly affecting key functional elements of miRNA precursors, such as DROSHA/DICER1 cleavage sites or the binding motifs of regulatory proteins. Depending on the effect of the variant, it may be considered a loss- or gain-of-function mutation. The further arguments that may support the potential functional/deleterious effect of mutations in miRNA genes are (i) the decreased frequency of SNPs (purifying selection) in these regions\u003csup\u003e\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e,\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e\u003c/sup\u003e, (ii) the results of testing artificial miRNA precursors/shRNAs, in which the tested sequence changes often resulted in altered processing of the precursors\u003csup\u003e\u003cspan additionalcitationids=\"CR30 CR31\" citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e–\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e\u003c/sup\u003e, and (iii) the results demonstrating the sequence specificity of DROSHA or DICER1 cleavage sites\u003csup\u003e\u003cspan additionalcitationids=\"CR34\" citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e–\u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e \u003cp\u003eTo shed more light on the consequences of mutations in miRNA genes (defined as sequences encoding pre-miRNAs with immediately adjacent 25 nt flanking sequences, roughly corresponding/overlapping with miRNA genes designated by the HUGO Gene Nomenclature Committee) on miRNA biogenesis, we took advantage of a large collection (n ~ 7000) of somatic mutations identified and annotated in miRNA genes\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e in cancer samples from The Cancer Genome Atlas (TCGA) project and the corresponding miRNA-seq data annotated with the use of isoMiRmap\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e. IsoMiRmap is a recently developed tool that allows precise annotation of miRNA reads considering length (isomiRs) and sequence (mutations) variation\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e. The comparison of the particular mutations with the corresponding miRNA-seq data allowed us to evaluate in real samples (not artificial functional models) the effect of the substantial number of mutations on miRNA levels, generated isomiRs, and the proportion of generated miRNA strands (5p/3p strand balance). Additionally, with the use of RNA-seq data for the selected mutations, we evaluated the effects of the mutations on the levels of the target mRNAs. The performed analyses revealed that a substantial number of tested mutations severely affect the function of miRNA genes and allow the identification of mutations with the most striking effects on particular aspects of biogenesis. This result is consistent with the idea that miRNA precursors are rather fragile structures in which even subtle changes affect their processing and suggests that most mutations in miRNA genes should be considered potentially functional (deleterious or gain of function) rather than neutral variants.\u003c/p\u003e \u003cp\u003eNotably, however, the mutations identified in this study as affecting miRNA genes should not be considered physiologically functional, nor do they necessarily play a role in cancer. Such analyses were beyond the scope of our study. By definition, only a small fraction of mutations occurring in the cancer genome (including mutations in miRNA genes) are expected to play a role in cancer and/or act as cancer-driving mutations. To prove such an effect, a specific set of functional tests focused on and adjusted to the expected impact of a particular mutation (or miRNA) in a specific cancer condition would have to be performed.\u003c/p\u003e "},{"header":"Methods","content":"\u003cp\u003e \u003cb\u003eData resources and generated metadata\u003c/b\u003e \u003c/p\u003e\u003cp\u003eThe IDs, sequences, and genomic coordinates of miRNAs (5p and 3p strands) and miRNA precursors were obtained from miRBase v.22.1\u003csup\u003e5\u003c/sup\u003e and used as references for mutation annotation and isomiRs classification. The miRNA gene IDs were used according to HUGO Gene Nomenclature. The following cancer type names and abbreviations were used according to TCGA nomenclature: adrenocortical carcinoma (ACC), bladder urothelial carcinoma (BLCA), breast invasive carcinoma (BRCA), cervical squamous cell carcinoma and endocervical adenocarcinoma (CESC), cholangiocarcinoma (CHOL), colon adenocarcinoma (COAD), lymphoid neoplasm diffuse large B-cell lymphoma (DLBC), esophageal carcinoma (ESCA), head and neck squamous cell carcinoma (HNSC), kidney chromophobe (KICH), kidney renal clear cell carcinoma (KIRC), kidney renal papillary cell carcinoma (KIRP), acute myeloid leukemia (LAML), brain lower grade glioma (LGG), liver hepatocellular carcinoma (LIHC), lung adenocarcinoma (LUAD), lung squamous cell carcinoma (LUSC), mesothelioma (MESO), ovarian serous cystadenocarcinoma (OV), pancreatic adenocarcinoma (PAAD), pheochromocytoma and paraganglioma (PCPG), prostate adenocarcinoma (PRAD), rectum adenocarcinoma (READ), sarcoma (SARC), skin cutaneous melanoma (SKCM), stomach adenocarcinoma (STAD), testicular germ cell tumor (TGCT), thyroid carcinoma (THCA), thymoma (THYM), uterine corpus endometrial carcinoma (UCEC), uterine carcinosarcoma (UCS), and uveal melanoma (UVM).\u003c/p\u003e\u003cp\u003eThe list of 7110 cancer somatic mutations in miRNA genes, together with sample identifiers and characteristics, was retrieved from Urbanek-Trzeciak et al. (Supplementary Table S2 in\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e). The mutations were identified in \u0026gt; 10,000 TCGA cancer samples representing 33 cancer types. To analyze the effect of mutations in miRNA genes on the generated miRNAs, we utilized small RNA-seq data of the corresponding TCGA samples mapped with the isoMiRmap tool\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e and retrieved from \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://cm.jefferson.edu/isoMiRmap/\u003c/span\u003e\u003cspan address=\"https://cm.jefferson.edu/isoMiRmap/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e \u003csup\u003e36\u003c/sup\u003e. IsoMiRmap allows precise annotation of miRNA length variants (isomiRs) and designated sequence variants, including miRNA reads with non-template nucleotides at the 3' end and miRNA reads with mutations (substitutions and 1-nt indels) annotated in Catalogue of Somatic Mutations in Cancer (COSMIC) v87 (released on 13 November 2018), covering a substantial fraction of the TCGA mutations\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e. For each sample, we downloaded three separate files generated by isoMiRmap, i.e., ‘exclusive’ (containing reads mapped to unique miRNA genes), ‘ambiguous’ (containing reads mapped to multiple genome sites/miRNA genes), and ‘snps’ (containing reads with genetic variants). A detailed description of isoMiRmap, the settings used for mapping, and the format of the output data can be found in\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e\u003cp\u003eTo prevent any uncertainties, biases, or confounding factors, the following mutations were removed from the analysis: (i) mutations identified in redundant miRNA genes (e.g., \u003cem\u003eMIR1-1\u003c/em\u003e (miR-1-1) and \u003cem\u003eMIR1-2\u003c/em\u003e); (ii) more than one mutation in a particular miRNA gene in one sample; (iii) long indels (more than 4 nt); (iv) mutations in samples with no small RNA-seq data; (v) mutations in miRNA genes with not defined one of the miRNA arms in miRBase; (vi) mutations in samples with more than one small RNA-seq dataset (sequenced multiple times); (vii) mutations in miRNA genes for which no reads were mapped in the mutated samples; (viii) mutations in samples with many ambiguous (mapping in multiple positions) reads in small RNA-seq; (ix) mutations located in mature miRNA or in proximity (± 2 nt) to DROSHA/DICER1 cleavage sites but not annotated in COSMIC v87.\u003c/p\u003e\u003cp\u003eOn the basis of the analysis of the above data and the above-listed criteria, we created a dataset consisting of 1309 mutations (1218 unique mutations) (Supplementary Table S1). Supplementary Table S1, in addition to data retrieved from\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e (yellowish columns), contains data generated in this study (white column). The new data were generated with the in-house script prepared in R. The list of mutations collected in Supplementary Table S1 served as a base list for selecting mutations for subsequent Experiments (see Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e1\u003c/span\u003e and Results).\u003c/p\u003e\u003cp\u003eTo analyze the effects of mutations on miRNA levels, isomiR distributions, strand balance, and the levels of miRNA targets (mRNAs) (analysis across multiple samples), we used crude small RNA-seq data retrieved from Loher et al.\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e and batch-corrected RNA-seq data retrieved from Hoadley et al.\u003csup\u003e\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e\u003c/sup\u003e, supplementary file \u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003eEBPlusPlusAdjustPANCAN_IlluminaHiSeq_RNASeqV2.geneExp.tsv\u003c/span\u003e. The small RNA-seq data were batch-corrected by us as previously described\u003csup\u003e\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e\u003c/sup\u003e, considering all isoMiRmap-annotated isomiRs. We retained all isomiRs to ensure that those potentially unique to mutated samples were not excluded. Batch correction was performed separately for each of the 33 TCGA cancer types (cohorts), considering the following confounding factors: (i) the platform (Illumina Genome Analyzer (GA) or Illumina HiSeq), (ii) tumor purity, (iii) and plate (in which the cDNA library was prepared). We did not correct the data for the type of RNA isolation protocol (direct (total RNA) or MultiMACs (poly-A depleted RNA)) because a single protocol was always used for a specific cancer type. Data about the platform and protocol were collected from the supplementary file in Hoadley et al.\u003csup\u003e\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e\u003c/sup\u003e and from the GDC Legacy Archive; tumor purity values were estimated using the TCGA tumor purity function from the TCGAbiolinks R package\u003csup\u003e\u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e\u003c/sup\u003e; and the plate IDs were identified on the basis of the aliquot barcodes. For batch correction, we used two algorithms: ComBat (version 3.80 from the sva R package\u003csup\u003e\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e\u003c/sup\u003e) or limma (version 3.6 from the limma R package\u003csup\u003e\u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e\u003c/sup\u003e available as part of the Bioconductor project). To identify potential confounders and select the proper algorithm, principal component analysis (PCA) was performed using the prcomp function of the stats R package (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.rdocumentation.org/packages/stats/versions/3.6.2/topics/prcomp.html\u003c/span\u003e\u003cspan address=\"https://www.rdocumentation.org/packages/stats/versions/3.6.2/topics/prcomp.html\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003e).\u003c/span\u003e On the basis of the PCA results, the following statistical tests were performed: the Kruskal‒Wallis test for categorical variables with n \u0026gt; 2, the Wilcoxon rank‒sum test for categorical variables with n = 2, and the Kendall tau rank correlation for ordinal variable purity. The effectiveness of the performed batch corrections was estimated using the FDR p-value, which indicates the association between potential confounders and the principal component (PC), and visualized with the PCA and t-distributed stochastic neighbor embedding (t-SNE) graphs (Rtsne wrapper for the Barnes–Hut t-Distributed Stochastic Neighbor Embedding from the Rtnse R package was used). On the basis of the above, limma was chosen for all tumor cohorts that required correction (CHOL, GBM, and UVM did not show any batch effects and thus were not batch-corrected).\u003c/p\u003e\u003cp\u003e \u003cb\u003emiRNA precursor structure modeling\u003c/b\u003e \u003c/p\u003e\u003cp\u003eThe secondary structures of the wild-type and mutant miRNA precursors and the change in Gibbs free energy (revised/optimized dG value) of the structures and individual base pairs or regions were predicted with the use of mfold version 3.6 or the mfold web server\u003csup\u003e\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e\u003c/sup\u003e with default parameters and processed in VARNA\u003csup\u003e\u003cspan citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e\u003c/sup\u003e. For modeling, we used pre-miRNA sequences (reconstructed based on miRBase as in\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e,\u003cspan citationid=\"CR44\" class=\"CitationRef\"\u003e44\u003c/span\u003e\u003c/sup\u003e) extended upstream (5’) and downstream (3’) by 25 nt flanking sequences. The structures with the lowest free energy are presented in this study. Spatial (3D) structures were modeled (on the basis of the predicted secondary structures) using RNAComposer software with the default parameters\u003csup\u003e\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e\u003c/sup\u003e and visualized via PyMOL\u003csup\u003e\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e\u003cp\u003e \u003cb\u003emiRNA level analysis\u003c/b\u003e \u003c/p\u003e\u003cp\u003eAnalysis of the influence of miRNA gene mutations on the miRNA level was performed using two different approaches. In the first approach [Experiment_1], within a specific mutant sample, we directly compared the levels (number of reads) of miRNA from the wild-type and the mutant alleles. To evaluate the fold change (FC) (or depletion) of the mutated allele, the numbers of reads at the RNA level were subsequently compared with the corresponding numbers of reads of mutant and wild-type alleles at the DNA level, and the significance of the difference was calculated with Fisher's exact test. This approach was applied only for mutations in mature miRNA sequences that can be distinguished at the RNA level. For the analysis, we considered only the ‘exclusive’ and ’mutated’ type of reads (wild-type and mutant). To select mutations adequate for the test and to avoid any uncertainties, the following mutations were excluded from the analysis: (i) not located in mature miRNA and (ii) located in miRNAs with fewer than 50 total miRNA reads in a mutated sample (not distinguishing reads from wild-type and mutant alleles). In the second approach [Experiment_3], in which we analyzed mutations located in all miRNA gene subregions (not only in mature miRNAs), we compared the level of miRNA in a mutated sample with the levels of miRNA in other samples of the same cancer type (without the mutation). For the analysis, we used batch-corrected RPM values. The level of a specific miRNA gene in a specific sample was calculated as the total RPM of all 5p and 3p isomiRs derived from a particular gene, including ‘exclusive’, ‘non-template’, and ‘mutated’ reads. To increase the normality of the miRNA level distribution, the RPM values were log2 transformed (RPM + 1). To avoid uncertainties, mutations with an RPM value lower than 10, either in a mutated sample or the average RPM of the other samples, were excluded. The level of a miRNA in a sample with a mutation vs. that in a in samples without mutations was compared with a one-sample t-test and z-score. To minimize the chance of false-positive results, only mutations with (i) a t-test with a Bonferroni-corrected p \u0026lt; 0.05, (ii) a z-score p \u0026lt; 0.05, (iii) |log2FC|\u0026gt;1, and (iv) identified as outliers defined as a datapoint \u0026gt; 1.5 IQR (interquartile range) below Q1 or above Q3 were considered significant.\u003c/p\u003e\u003cp\u003e \u003cb\u003eIsomiR classes\u003c/b\u003e \u003c/p\u003e\u003cp\u003eConsistently with the convention described previously\u003csup\u003e\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e\u003c/sup\u003e, all isomiRs were categorized into 9 basic classes denoted as follows: 0|0 (reference (canonical) isomiR as annotated in miRBase); 0|-; 0|+; +|0; -|0; +|-; -|+; +|+; -|-. The signs/values before and after the vertical line correspond to the 5p and 3p miRNA ends (5p|3p) and indicate the direction of the particular end shift, i.e., no change (0), upstream (3p to 5p, (-)), and downstream (+). Additionally, when comparing isomiR profiles between samples (mutant vs. wild-type), reads with non-template 3p-end modifications were considered an additional class of (non-template) isomiRs, denoted (nt) if any modification was present. The extension of the classification may be used to indicate the exact isomiR coordinates, e.g., -2|-1 (indicating an isomiR with the 5p-end shifted by 2 nucleotides upstream and the 3p-end shifted by 1 nucleotide upstream) or by using (n) for any nucleotide change at the 3p-end to emphasize the modifications at the 5p-end. The principles of isomiR classification and denotation are shown in Supplementary Figure S1.\u003c/p\u003e\u003cp\u003e \u003cb\u003eIsomiR distribution analysis\u003c/b\u003e \u003c/p\u003e\u003cp\u003eAs in the case of the analysis of the effect of mutations on the level of miRNA, the effect of mutations on isomiR distribution was analyzed using two different approaches. In the first approach, we compared the fraction of isomiRs (isomiR profiles) of mutant and wild-type alleles, which were divided into 9 isomiR classes (as defined above), within one sample. To calculate the number and fraction of reads classified into particular isomiR classes, we used raw read counts (only the ‘exclusive’ and “mutated” type of reads were taken into account). Rare cases of isomiRs/reads for which it could not be determined whether they originated from wild-type or mutant alleles were excluded from the analysis. From the list of mutations (Supplementary Table S1), we excluded mutations that (i) are not located in mature miRNAs and (ii) have \u0026lt; 50 reads derived from either allele. For statistical analysis of changes in isomiR profiles, we used Pearson’s chi-square independence test and Cramer’s V test. For the chi-square test, we added a pseudocount (number of reads + 5) to each isomiR class. Cramer’s V value measures the strength of the association/relationship between two variables. Cramer’s V values \u0026lt; 0.1, 0.1–0.2, 0.2–0.4, and \u0026gt; 0.4 were interpreted as no association or small (weak), medium, and large (strong) associations, respectively\u003csup\u003e\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e\u003cp\u003eIn the second approach, we compared the isomiRs (isomiR profiles) derived from a mutated gene in the sample with the mutation and in other samples of the same cancer type without the mutation. From the list of mutations (Supplementary Table S1), we excluded mutations that (i) had less than 20 RPM for a particular miRNA strand in the mutated sample and less than 20 RPM for the average in the corresponding wild-type samples. If the criteria were fulfilled for both strands, the isomiR profiles were analyzed for each strand. As described in the subchapter above, isomiRs were divided into 10 classes, and fractions of particular isomiR classes were calculated on the basis of batch-corrected RPM values. For each miRNA of interest, we subsequently calculated the distance (sum of differences of all isomiR fractions; values from 0 to 2 where 0 indicates no difference and 2 indicates the occurrence of completely different isomiRs) between the isomiR profile in a sample with mutation and an average isomiR profile of samples without mutation. The distance calculated for the sample with the mutation was compared (one-sample t-test and z-score) with similarly calculated distances of all samples without the mutation and visualized on a cumulative graph of isomiR profiles. Only mutations with (i) a t-test with a Bonferroni-corrected p \u0026lt; 0.05, (ii) a z-score p \u0026lt; 0.05, and (iii) a distance from the average isomiR profile of samples without a mutation \u0026gt; 0.2 were considered significant.\u003c/p\u003e\u003cp\u003e \u003cb\u003eStrand balance analysis\u003c/b\u003e \u003c/p\u003e\u003cp\u003eStrand balance analysis was performed by comparing the 5p/3p strand balance (calculated as log2(5p_raw_counts + 1/3p_raw_counts + 1), including ‘exclusive’, ‘non-template’, and ‘mutated’ reads) of a mutated miRNA gene in a mutated sample vs. the strand balance of the gene in other samples of the same cancer type without mutation. To select mutations that were adequate for the test and to avoid any uncertainties, mutations with fewer than 50 total reads for a particular miRNA gene (Supplementary Table S1) were excluded from the analysis. To ensure a reliable estimation of the strand balance in the other (reference) samples, only samples with a total number of reads ≥ 50 were included in the analysis. To identify mutations that significantly affect strand balance, the values were compared with a one-sample t-test and z-score. Only mutations with (i) a t-test with a Bonferroni-corrected p \u0026lt; 0.05, (ii) a z-score p \u0026lt; 0.05, (iii) \u0026gt; 10% change in strand fractions between mutated sample and average from other WT samples, and (iv) identified as outliers, defined as a datapoint \u0026gt; 1.5 IQR below Q1 or above Q3, were considered significant.\u003c/p\u003e\u003cp\u003e \u003cb\u003eDifferential expression\u003c/b\u003e \u003c/p\u003e\u003cp\u003eFor differential expression analysis, we used batch-corrected mRNA (RNA-seq) data (see above in section ‘Data resources and generated metadata’). The analysis was performed with the use of the DESeq2-based\u003csup\u003e\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e\u003c/sup\u003e R pipeline (SARTools)\u003csup\u003e\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e\u003c/sup\u003e to compare selected samples with mutations affecting the analyzed aspects of miRNA gene biogenesis with other samples of the same cancer type without a mutation.\u003c/p\u003e\u003cp\u003e \u003cb\u003emiRNA target prediction\u003c/b\u003e \u003c/p\u003e\u003cp\u003emiRNA target prediction was performed using TargetScan\u003csup\u003e\u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e\u003c/sup\u003e and mirDB\u003csup\u003e\u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e\u003c/sup\u003e, allowing target prediction for both wild-type and mutant (custom) sequences or isomiRs. For the prediction of mutation-specific targets, we used the most common and/or the most differentiated isomiRs in a mutated sample. We subsequently compared the list of predicted targets with the list of differentially expressed genes between the wild-type and mutant samples.\u003c/p\u003e\u003ch2\u003eStatistical analysis\u003c/h2\u003e\u003cp\u003eStatistical analyses were performed using the stats R package or the Real Statistics Resource Pack for MS Excel (\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://real-statistics.com/\u003c/span\u003e\u003cspan address=\"https://real-statistics.com/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003cspan type=\"Underline\" class=\"Underline\" name=\"Emphasis\"\u003e).\u003c/span\u003e\u003c/p\u003e"},{"header":"Results","content":"\u003cp\u003e \u003cb\u003eData processing\u003c/b\u003e \u003c/p\u003e \u003cp\u003eTo evaluate the effects of mutations in miRNA genes on the levels and sequences of generated miRNAs, we compiled data from the list of 7,110 cancer somatic mutations identified in miRNA genes in \u0026gt;\u0026thinsp;10,000 TCGA cancer samples representing 33 cancer types\u003csup\u003e\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e\u003c/sup\u003e and from the corresponding miRNA reads (small RNA-seq datasets) mapped using isoMiRmap\u003csup\u003e\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e\u003c/sup\u003e. IsoMiRmap is a tool that, in addition to mapping wild-type miRNAs, allows reads mapping to mutations identified in TCGA cancer samples, annotated in Catalogue of Somatic Mutations in Cancer (COSMIC) v87\u003csup\u003e36\u003c/sup\u003e (for details, see Methods). To avoid confusion resulting from the ambiguous mapping of miRNA reads, we excluded mutations located in redundant/duplicated miRNA genes [e.g., \u003cem\u003eMIR1-1\u003c/em\u003e (miR-1-1) and \u003cem\u003eMIR1-2\u003c/em\u003e or \u003cem\u003eMIRLET7A1\u003c/em\u003e (let-7a-1), \u003cem\u003eMIRLET7A2\u003c/em\u003e, and \u003cem\u003eMIRLET7A3\u003c/em\u003e] from the analysis. Additionally, we removed mutations/samples for which there was no small RNA-seq data and mutations in genes for which no miRNA reads were detected in the corresponding samples. Other exclusion criteria are described in the Methods section (\u0026ldquo;Data resources and generated metadata\u0026rdquo;) and presented in Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e1\u003c/span\u003e. The remaining 1309 mutations listed/characterized in Supplementary Table S1 served as a resource of mutations for subsequent analyses. The effect of mutations on miRNAs expressed from the mutated genes was analyzed using two general approaches, i.e., comparing miRNAs expressed from a particular miRNA gene from the mutant and wild-type alleles in the mutated sample (APPROACH_1) and comparing miRNAs expressed from a mutated gene in the sample with the mutation with other samples of the same cancer type without the mutation (APPROACH_2). The general concept of the study, as well as the exclusion criteria and the number of mutations selected for subsequent Experiments, are shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig7\" class=\"InternalRef\"\u003e1\u003c/span\u003e and Supplementary Figure S2.\u003c/p\u003e \u003cp\u003e \u003cb\u003eAnalysis of the effect of miRNA gene mutations by direct comparison of miRNAs generated from mutated and wild-type alleles in individual samples with mutations\u003c/b\u003e \u003c/p\u003e \u003cp\u003eIn this approach, we took advantage of the fact that some mutations may be directly observed at the RNA level; therefore, direct observation of the miRNAs derived from the mutant and wild-type alleles in one sample is possible without biases resulting from the biological and technical differences between samples or their processing/analysis. The limitation of such an approach, however, is that it allows the analysis only of mutations located in sequences encoding mature miRNAs (expressed at the miRNA level). In this approach, to ensure proper distinction between wild-type and mutant reads, only mutations annotated in COSMIC v87 were taken into account.\u003c/p\u003e \u003cp\u003eIn the first experiment (Experiment_1), to investigate the effects of mutations on miRNA levels, we compared the allelic fractions of the mutations (proportions of reads with mutations) at the genomic (DNA) and transcript (mature miRNA) levels. To ensure proper determination of the mutant fraction, only mutations in miRNAs with at least 50 miRNA reads (wild-type\u0026thinsp;+\u0026thinsp;mutant) were used for analysis (n\u0026thinsp;=\u0026thinsp;53). Surprisingly, as shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e2\u003c/span\u003eA, there was very little relationship between the fractions of mutant alleles at the DNA and RNA levels, and most mutations (n\u0026thinsp;=\u0026thinsp;32, 60%) significantly (adjusted p\u0026thinsp;\u0026lt;\u0026thinsp;0.05; Fisher\u0026rsquo;s exact test) deviated from the trendline (x\u0026thinsp;=\u0026thinsp;y), representing an equal proportion of mutant reads at the DNA and miRNA levels, expected for mutations with no effect on the miRNA level. Only one mutation significantly increased the miRNA level (adj. p\u0026thinsp;\u0026gt;\u0026thinsp;0.05; absolute FC\u0026thinsp;\u0026ge;\u0026thinsp;2), whereas the vast majority of these mutations decreased the miRNA level of mutated alleles (1 vs. 31; p\u0026thinsp;=\u0026thinsp;0.00003; Fisher\u0026rsquo;s exact test), in most cases almost to 0. Moreover, the relative decreases (FCs) had much larger amplitudes than the relative increases in miRNA levels (Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e2\u003c/span\u003eA inset). As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e2\u003c/span\u003eB, the mutations affecting miRNA levels are roughly equally distributed along mature miRNA sequences and do not cluster in any specific position of the precursor (details regarding each of the mutations are listed in Supplementary Table S2).\u003c/p\u003e \u003cp\u003eIt is noteworthy that several mutations (n\u0026thinsp;=\u0026thinsp;5) that occurred in more than one sample, especially in samples of the same cancer type, generally had consistent effects on miRNA levels (Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e2\u003c/span\u003eA, inset), further confirming the validity of the analysis. Among these mutations, n.22G\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR379\u003c/em\u003e was identified in two LUAD samples and, in both cases, acutely decreased the level of mutant miRNA, whereas n.75C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR342\u003c/em\u003e (in 2 SKCM samples), n.59T\u0026thinsp;\u0026gt;\u0026thinsp;C in \u003cem\u003eMIR142\u003c/em\u003e (in 2 DLBC samples), and n.51G\u0026thinsp;\u0026gt;\u0026thinsp;A in \u003cem\u003eMIR320A\u003c/em\u003e (in LUAD, LUSC, and HNSC samples) had no significant effect on the miRNA level. Only n.35C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR205\u003c/em\u003e had a discordant effect on the miRNA level in the SKCM sample (decrease) vs. the CESC and LUSC samples (no significant change). These discrepancies may, however, result from differences in cancer/tissue type, differences in mutation allelic frequency, and general variation in cancer samples.\u003c/p\u003e \u003cp\u003eRecently, it was shown that the nucleotide composition of a 3-base pair-long motif of the miRNA precursor duplex encompassing the DICER1 cleavage site, called GYM, affects the efficiency of DICER1 processing and thus may affect miRNA levels\u003csup\u003e\u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e53\u003c/span\u003e\u003c/sup\u003e. By comparing the GYM scores of the 3 mutations located in the GYM region with those of their WT counterparts, we found that the effect of one of the mutations, i.e., n.35C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR523\u003c/em\u003e, which decreases the miRNA level, may be well explained by the GYM effect, i.e., the GYM score decreased from 26 for the WT precursor to 10 for the mutant precursor (the GYM score ranges from 0-100, with higher values indicating more efficient processing; Fig.\u0026nbsp;\u003cspan refid=\"Fig8\" class=\"InternalRef\"\u003e2\u003c/span\u003eC). Notably, however, miRNA levels can be influenced by many other factors, including the structure of the miRNA precursor (discussed below).\u003c/p\u003e \u003cp\u003eIn Experiment_2, to investigate the impact of mutations on the precision of DROSHA/DICER1 cleavage, i.e., generated isomiRs, we compared the isomiR profiles of particular miRNAs generated from wild-type and mutant alleles. To ensure the reliable determination of isomiR distributions, for the analysis, we selected only mutations/samples with at least 50 miRNA reads for both the wild-type and mutant alleles (n\u0026thinsp;=\u0026thinsp;16; Supplementary Table S3). To compare the distribution of the corresponding wild-type and mutant miRNAs, we classified all the miRNA reads into 9 isomiR classes categorized on the basis of the position of their ends (5p|3p) upstream (+) or downstream (-), against the corresponding ends annotated in miRBase (0|0) (as proposed previously\u003csup\u003e\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e\u003c/sup\u003e and graphically illustrated in Supplementary Figure S1). Ten of the 16 (63%) tested mutations induced significant (chi-square adj. p\u0026thinsp;\u0026lt;\u0026thinsp;1E-11) changes in isomiR profiles and met the criterion of Cramer\u0026rsquo;s V\u0026thinsp;\u0026gt;\u0026thinsp;0.2, indicating at least a moderate relationship of the mutations with the isomiR profile (Supplementary Table S3). Among the identified mutations, 3 indicate a strong effect on the isomiR profile (Cramer\u0026rsquo;s V\u0026thinsp;\u0026gt;\u0026thinsp;0.4). The mutations affecting isomiR profiles are roughly equally distributed along mature miRNA sequences (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eA). To determine whether Cramer\u0026rsquo;s V value reliably distinguishes changes induced by mutations from random variation in isomiR distributions, we also compared isomiR profiles of the same mutants and wild-type alleles of the same miRNAs in different samples. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eB, the Cramer\u0026rsquo;s V values are much lower for pairs of the same mutants or wild-type alleles than for pairs of corresponding mutant and wild-type alleles from the same samples, even though the former come from different samples or even different cancer types. Thus, a Cramer's V threshold of 0.2 accurately distinguishes significant changes from random variability. Most of the identified mutations induce changes of \u0026gt;\u0026thinsp;2-fold in at least one isomiR class and, in most cases, substantially (\u0026gt;\u0026thinsp;20%) change the level of the main isomiRs (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eC, D, E; Supplementary Figure S3). Notably, 6 of the mutations induced changes predominantly at the 5p-end, resulting in a shift in the miRNA seed sequences, which have a direct impact on target recognition.\u003c/p\u003e \u003cp\u003eAmong the mutations with the strongest effect on the isomiR distribution was n.75C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR342\u003c/em\u003e, located at the 15th nucleotide of miR-342-3p, which shortened the miRNA at the 3p-end and thus severely reduced (\u0026gt;\u0026thinsp;4-fold) the fraction of the canonical (0|0) and 0|+ isomiRs in favor of the 0|- isomiR (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eC). This mutation was identified in two SKCM samples. The very similar isomiR profile in two independent samples with this mutation (Cramer\u0026rsquo;s V\u0026thinsp;\u0026gt;\u0026thinsp;0.2) proves the validity of the results and the specificity of the changes induced by the mutation.\u003c/p\u003e \u003cp\u003eAn example of a mutation affecting the 5p-end of the generated miRNAs is n.55A\u0026thinsp;\u0026gt;\u0026thinsp;G in \u003cem\u003eMIR142\u003c/em\u003e, located at the 3rd nucleotide of the miR-142-3p seed (Cramer\u0026rsquo;s V\u0026thinsp;=\u0026thinsp;0.638). This mutation elongates the 5p-end of the generated isomiRs by one nucleotide, severely reducing the level of the +|- and +|0 isomiRs (predominant in the wild-type allele) and increasing the level of the 0|- and 0|0 isomiRs (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eD). A similar shift (in favor of isomiRs that are longer at the 5p end) is induced by 3 other mutations located at nearby positions in the miR-142-3p seed (n.58G\u0026thinsp;\u0026gt;\u0026thinsp;C and n.59T\u0026thinsp;\u0026gt;\u0026thinsp;C (2x; in 2 samples), although the strength of the effect varies depending on the mutation and is not formally significant (Cramer\u0026rsquo;s V\u0026thinsp;~\u0026thinsp;0.18) in all cases (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eE). We subsequently combined all the isomiRs of the \u003cem\u003eMIR142\u003c/em\u003e mutant samples on the basis of particular changes at the 5p-end (5p-isomiRs; ignoring changes at the 3p-end). As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eF, the predominant isomiR class expressed from the wild-type allele (accounting for approximately 80% of all 4 samples) consists of isomiRs\u0026thinsp;+\u0026thinsp;1|n shortened by 1 nt at the 5p-end. The +\u0026thinsp;1|n isomiRs are expressed at much lower levels in the mutant alleles, ranging from 9 to 61%, depending on the mutation. In contrast, the mutant isomiRs consisted of a greater fraction (38\u0026ndash;89%) of isomiRs with a canonical 5p-end (0|n), accounting for only\u0026thinsp;~\u0026thinsp;20% of the wild-type isomiRs (Fig.\u0026nbsp;\u003cspan refid=\"Fig9\" class=\"InternalRef\"\u003e3\u003c/span\u003eF). This shift of the 5p-miRNA end (+\u0026thinsp;1|n to 0|n) results in an additional change in the seed sequence that is independent of the point changes directly introduced by the individual mutations. The remaining mutations that significantly affect isomiR distribution are shown in Supplementary Figure S3.\u003c/p\u003e \u003cp\u003e \u003cb\u003eAnalysis of the effect of miRNA gene mutations by comparing miRNAs in samples with the mutation versus corresponding samples without the mutation\u003c/b\u003e \u003c/p\u003e \u003cp\u003eTo extend the analysis to mutations located in other parts of miRNA genes (not only in mature miRNAs that can be observed at the RNA level), we used another approach, in which we analyzed the effect of mutations by comparing miRNAs expressed from a mutated gene in a sample with the mutation with other samples of the same cancer type that lack the mutation. However, as the power of such an approach is strongly limited by the dilution of the mutation effect by the presence of the normal allele and the contamination of cancer samples with normal (non-cancerous) cells (in most cases, a mutant allele accounts for \u0026lt;\u0026thinsp;\u0026lt;\u0026thinsp;50%; Supplementary Table S1), as well as very high genetic and transcriptional variation in cancer samples, we sacrificed sensitivity to focus on identifying single mutations with the most profound effects. For this purpose, we used the above standard multiple criteria for classifying results as significant (see Methods), minimizing the potential of false-positive results but likely missing the effects of many mutations (false-negatives).\u003c/p\u003e \u003cp\u003eIn Experiment_3, to identify mutations that affect the level of miRNAs, we compared batch-corrected levels of miRNAs (log2-transformed RPM values of the 5p and 3p arms) in samples with mutations to the average levels of the miRNAs in samples without mutations. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003eA, of the 682 mutations (648 unique mutations) that fulfilled the criteria of the analysis, 21 significantly affected the miRNA level, including 15 mutations that increased the level and 6 mutations that decreased the level (Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003e, Supplementary Table S4). The observed excess of mutations increasing the level of miRNA (compared to Experiment_1) is likely attributable to the reduced statistical power to detect downregulation relative to upregulation. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003eB, mutations affecting miRNA levels are located in all subregions of the miRNA precursors and are roughly equally distributed along the sequence. Among the examples of mutations affecting miRNA levels are mutation n.5-7delAGC in the 5p-flank in \u003cem\u003eMIR122\u003c/em\u003e, which causes a decrease in the miRNA level (log2(FC)\u0026lt;-6) (Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003eC), and mutation n.85G\u0026thinsp;\u0026gt;\u0026thinsp;A, which is located in the 3p-flank in \u003cem\u003eMIR518E\u003c/em\u003e, which leads to a dramatic increase in the miRNA level (log2(FC)\u0026thinsp;\u0026gt;\u0026thinsp;10) (Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003eD). Both mutations affect the predicted 2D structure and stability of the miRNA precursors by relaxing the duplex structure near the mutation site. Moreover, both mutations alter the expression of numerous genes, including the predicted targets of the relevant miRNAs (Fig.\u0026nbsp;\u003cspan refid=\"Fig10\" class=\"InternalRef\"\u003e4\u003c/span\u003eE, F; Supplementary Table S5).\u003c/p\u003e \u003cp\u003eIn Experiment_4, to investigate the effect of mutations on the precision of DROSHA/DICER1 cleavage, for each mutation, we compared the isomiR profiles of a given miRNA in the sample with the mutation and corresponding (same cancer type) samples without mutations. To ensure the reliable determination of isomiR profiles, only samples with at least 20 RPM were considered in the analysis. Among the 420 mutations (397 unique mutations) that fulfilled the criteria, 32 mutations significantly affected the isomiR profile in at least 1 strand (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e5\u003c/span\u003eA, Supplementary Table S6). Mutations affecting the isomiR profile are located in all subregions of the miRNA precursors and are roughly equally distributed along the sequence (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e5\u003c/span\u003eB).\u003c/p\u003e \u003cp\u003eAmong the mutations affecting the isomiR profile is n.85G\u0026thinsp;\u0026gt;\u0026thinsp;A, which is located in the 3p-flank of \u003cem\u003eMIR142\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e5\u003c/span\u003eC). The mutation leads to extension by one nucleotide of the 5p-end of the generated miR-142-3p isomiRs, severely reducing the level of the +|0 and +|- isomiRs (predominant in the wild-type samples) and increasing the level of the 0|0 and 0|- isomiRs. Interestingly, other mutations located in different parts of \u003cem\u003eMIR142\u003c/em\u003e induced a similar shift in the isomiR profile (as analyzed in Experiment_2, compare the isomiR profile graph in Figure ). Another mutation that alters the isomiR profile is n.21-22delTA in \u003cem\u003eMIR10B\u003c/em\u003e, which is located in the 5p flank of the gene (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e5\u003c/span\u003eD). The mutation causes the extension of the 5p-end of miR-10b-5p by one or two nucleotides, leading to a reduction in the fraction of 0|n isomiRs in favor of +\u0026thinsp;1|n and +\u0026thinsp;2|n isomiRs. Another example is n.-6C\u0026thinsp;\u0026gt;\u0026thinsp;G in \u003cem\u003eMIR21\u003c/em\u003e, which is located in the 5p flank of the gene (Fig.\u0026nbsp;\u003cspan refid=\"Fig11\" class=\"InternalRef\"\u003e5\u003c/span\u003eE). The mutation induces shortening by 1 nt of the 3p miRNA at its 3p-end, severely reducing the fraction of 0|+ isomiRs predominant in WT samples in favor of the canonical (0|0) isomiR. As shown in the volcano plots (Supplementary Figure S4), some genes are differentially expressed in the mutated samples. However, the relationship between the differentially expressed genes and the observed isomiR shifts cannot be directly determined. The isomiR profiles of other mutations that significantly affect isomiR profiles are presented in Supplementary Figure S5.\u003c/p\u003e \u003cp\u003eIn Experiment_5, we analyzed the effect of mutations on miRNA strand balance, i.e., the ratio of miRNAs derived from 5p and 3p miRNA arms of precursor, compared between mutant and wild-type samples as log2(5p\u0026thinsp;+\u0026thinsp;1/3p\u0026thinsp;+\u0026thinsp;1). To ensure reliable determination of the strand balance, only samples with \u0026ge;\u0026thinsp;50 reads of a given miRNA were considered. As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eA, from the 515 mutations (490 unique mutations) fulfilling the criteria of the analysis (for details, see Methods), 4 mutations significantly affected the miRNA strand balance, including one mutation located in each of the following regions: 5p-flank, 3p-flank, loop, and miRNA duplex (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eA, B; Supplementary Table S7).\u003c/p\u003e \u003cp\u003eThe most striking example is the duplication of two Gs (n.97_98dupGG) in the 3p-flank in \u003cem\u003eMIR205\u003c/em\u003e identified in THYM. In all the THYM samples, miR-205-5p was the predominant strand, accounting for more than 99.7% of all the reads (similar to other TCGA cancers, as well as other tissues reported in different databases). In contrast, in the mutant sample, the fraction of miR-205-5p decreased substantially (down to 65%), resulting in a dramatic increase (\u0026gt;\u0026thinsp;100x, up to 35%) of the miR-205-3p fraction (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eC). It is noteworthy that the observed change of the 5p/3p balance is detectable even though the fraction of the mutant allele (at the DNA level) accounts for only 19%, which suggests that the actual effect of the mutant allele is much stronger (this comment also applies to other mutations mentioned in this section). A closer look at the miRNAs generated from the 5p and 3p arms revealed that, in addition to changing the 5p/3p balance, the mutation also shifted (by 1 to 3 nucleotides toward the base of the precursor hairpin) the predominant DROSHA/DICER1 cleavage sites, generating abnormal miRNAs (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eC). To determine whether this dramatic strand/isoform imbalance (the mutant-specific miRNA) affects gene expression, we performed differential expression analysis and identified 31 and 5 genes whose expression was downregulated and upregulated, respectively, in the sample with the \u003cem\u003eMIR205\u003c/em\u003e mutation. Among the downregulated genes, four targets (predicted by TargetScan and/or mirDB) of mutant-specific miR-205-5p (-2|-1) and 3p (+\u0026thinsp;2|+3) were identified (Supplementary Table S5; Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eE).\u003c/p\u003e \u003cp\u003eAn example of a mutation that disturbs the balance in the opposite direction (in favor of the 5p strand) is n.84C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR365B\u003c/em\u003e, which is located in the 3p-strand of the miRNA duplex (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eD). This mutation increases the level of the 5p-strand more than fivefold, from 3\u0026ndash;19%. As the mutation destabilizes the miRNA duplex at the 5p-end of miR-365b-5p (replaces G:C with the G:U pair and increases the duplex free energy by 2 kcal/mol), its effect of favoring the 5p over the 3p strand is consistent with the principle of greater preference/stability (more efficient loading into the miRISC complex) of the strand with a more relaxed (less stable) 5p-end\u003csup\u003e\u003cspan citationid=\"CR54\" class=\"CitationRef\"\u003e54\u003c/span\u003e\u003c/sup\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eD). As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig12\" class=\"InternalRef\"\u003e6\u003c/span\u003eF, many genes were downregulated in the mutated sample, and 7 of these genes were predicted targets of miR-365b-5p, the fraction of which was increased in the mutant sample. The other two mutations that affect the strand balance are shown in Supplementary Figure S6.\u003c/p\u003e \u003cp\u003e \u003cb\u003eThe impact of miRNA gene mutations on the structure of miRNA precursors\u003c/b\u003e \u003c/p\u003e \u003cp\u003eAlthough mutations can affect the functioning of miRNA genes through various mechanisms, the most direct/common impact seems to be modifying the structure of miRNA precursors. Therefore, to preliminarily estimate the effects of mutations on the structure of miRNA precursors, we modeled and compared the wild-type and mutant precursor structures of all the tested miRNA genes (n\u0026thinsp;=\u0026thinsp;1309). Most mutations increase the change in the Gibbs free energy (dG) value between mutant and WT precursors, i.e., decrease the stability of the precursor structures (Supplementary Table S8). The average dG value of the mutant structures was 0.77 kcal/mol greater than that of the corresponding wild-type structures (p\u0026thinsp;=\u0026thinsp;1.5E-35; paired t-test), and there was an excess of mutations increasing the dG value (destabilizing the precursor structure) over mutations decreasing the dG value (Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e7\u003c/span\u003eA). Further comparisons revealed that the increase in dG (ddG) was significantly greater for mutations affecting the functioning of miRNA genes (functional mutations identified in Experiments_1\u0026ndash;5) than for the remaining mutations (t-test, p\u0026thinsp;=\u0026thinsp;1.73E-08). As shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e7\u003c/span\u003eB, the functional mutations cluster at the top of the ddG ranking. The ddG for functional mutations is independent of the functional effect of the mutation; it applies both to mutations impacting miRNA levels (t-test, p\u0026thinsp;=\u0026thinsp;6.92E-09) and mutations affecting isomiR distributions (t-test, p\u0026thinsp;=\u0026thinsp;0.0012) (Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e7\u003c/span\u003eB; due to the low number of mutations, we did not test mutations affecting 5p/3p-strand balance). The relationship between the results of structural and functional analyses further confirms the reliability of the identified functional mutations (results of Experiments_1\u0026ndash;5), as both analyses are completely independent. Examples of mutations affecting the miRNA precursor structures are shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig13\" class=\"InternalRef\"\u003e7\u003c/span\u003eC.\u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003eThanks to the genetic code, predicting the consequences of mutations and distinguishing deleterious from neutral mutations in coding sequences is relatively straightforward (although not trivial). On the other hand, owing to modern sequencing approaches, such as WGS, thousands of new genetic variants, primarily located in noncoding sequences, are being identified in each sequenced sample. Unfortunately, due to the heterogeneity of noncoding elements and a lack of appropriate \u0026quot;codes\u0026quot; (rules), predicting the consequences of mutations or even estimating the likelihood of deleteriousness for mutations in the best-defined noncoding elements, such as miRNA genes, is nearly impossible. Consequently, the vast majority of sequence variants identified in the noncoding genome are assumed to be neutral or variants of unknown significance (VUS) and are ignored. To overcome these limitations and ultimately develop \u0026quot;codes\u0026quot; to predict the consequences of variants in noncoding elements, more data need to be collected on the impacts of mutations on the functionality of these elements.\u003c/p\u003e\n\u003cp\u003eAmong the most well-recognized and highly conserved noncoding elements are miRNA genes, particularly their ~\u0026thinsp;100 nt long regions encoding hairpin-structured miRNA precursors (pre-miRNAs). Although analysis of genetic variation in the noncoding genome has long been neglected, numerous variants in miRNA genes have been extracted from whole-genome datasets. Moreover, few miRNA gene mutations have been implicated in Mendelian diseases, including mutations in \u003cem\u003eMIR96\u003c/em\u003e in non-syndromic hearing loss\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e11\u003c/span\u003e\u003c/sup\u003e, mutations in \u003cem\u003eMIR184\u003c/em\u003e in different hereditary eye diseases\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e12\u003c/span\u003e\u003c/sup\u003e, mutations in \u003cem\u003eMIR204\u003c/em\u003e in retinal dystrophy\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e13\u003c/span\u003e\u003c/sup\u003e, and mutations in \u003cem\u003eMIR140\u003c/em\u003e in skeletal dysplasia\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e14\u003c/span\u003e\u003c/sup\u003e. There are also a few miRNA genes whose mutations have been suggested to be potential cancer drivers, including the \u003cem\u003eMIR15A/MIR16-1\u003c/em\u003e cluster in CLL, \u003cem\u003eMIR142\u003c/em\u003e in different blood cancers, \u003cem\u003eMIR122\u003c/em\u003e in liver carcinoma, and \u003cem\u003eMIR21\u003c/em\u003e in various cancers\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e7\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e15\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e55\u003c/span\u003e\u0026ndash;\u003cspan class=\"CitationRef\"\u003e58\u003c/span\u003e\u003c/sup\u003e. Nevertheless, a larger-scale analysis of the molecular consequences of miRNA gene variants has never been performed, and almost all previous studies focused on mutations identified in the seed sequences\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e11\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e13\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e14\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e17\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e18\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e59\u003c/span\u003e\u0026ndash;\u003cspan class=\"CitationRef\"\u003e62\u003c/span\u003e\u003c/sup\u003e. Only a few studies have investigated the impact of variants in other regions of miRNAs\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e19\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e20\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e25\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e63\u003c/span\u003e\u0026ndash;\u003cspan class=\"CitationRef\"\u003e67\u003c/span\u003e\u003c/sup\u003e (summarized in\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e9\u003c/span\u003e\u003c/sup\u003e). Moreover, almost all of these studies have investigated only the impact on the miRNA level. To date, the only and largest systematic analysis of sequence variants in miRNA genes, performed more than 15 years ago, has examined the effects of 24 common SNPs located in different parts of miRNA precursors. Using several cellular functional assays, the study revealed that most of the variants disturb the proper function of miRNA genes, including the efficiency of miRNA generation and precision of miRNA processing\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e20\u003c/span\u003e\u003c/sup\u003e.\u003c/p\u003e\n\u003cp\u003eTherefore, taking advantage of the large collection of mutations and corresponding transcriptome (miRNA-seq) datasets generated by the TCGA project, we performed the most comprehensive analysis to date of the impact of mutations on miRNA gene functioning, exceeding the number of variants studied in previous projects by orders of magnitude. In contrast to most previous studies, the approach allowed us, for the first time, to analyze the effects of mutations in their real genetic context, under natural conditions in real cells/tissues, without artificially generated models using engineered genes, overexpression, or generic cell lines. Furthermore, in contrast to most previous studies, our analysis included mutations located in all parts of the miRNA gene, not only in the seed region. An additional advantage of our study is that the vast majority of somatic mutations are randomly occurring variants. This ensures that the analysis is not biased toward more likely neutral variants, such as common SNPs, or deleterious variants, such as mutations identified in Mendelian diseases. Moreover, in addition to the analysis of the impact of mutations on the miRNA level, which has been mostly analyzed in other studies, we also analyzed the effects of mutations on the miRNA strand balance and isomiR profiles.\u003c/p\u003e\n\u003cp\u003eAltogether, we analyzed the effects of 703 mutations selected in multiple steps from ~\u0026thinsp;7000 mutations identified in miRNA genes in TCGA samples (Fig. \u003cspan class=\"InternalRef\"\u003e1\u003c/span\u003e, Supplementary Figure S2). The mutations were selected on the basis of very strict criteria, ensuring that their effects may be reliably evaluated (e.g., by excluding mutations in redundant miRNA genes such as \u003cem\u003eMIR1-1\u003c/em\u003e and \u003cem\u003eMIR1-2\u003c/em\u003e). Nevertheless, this is the largest and most complex analysis of miRNA gene mutations performed to date.\u003c/p\u003e\n\u003cp\u003eFirst, by comparing the levels of mutant and wild-type alleles at the genomic (DNA) and transcriptomic (miRNA) levels, we showed that most miRNA gene mutations (at least those located in the miRNA duplex) affect the miRNA level (32/53) and/or precision of DROSHA/DICER1 cleavages (isomiR profile 10/16), i.e., they are deleterious for the proper functioning of miRNA genes at the molecular level. Additional mutations, which are located in the whole pre-miRNA and affect miRNA gene functionalities, including miRNA levels (n\u0026thinsp;=\u0026thinsp;21), isomiR profiles (n\u0026thinsp;=\u0026thinsp;32), and miRNA strand balance (n\u0026thinsp;=\u0026thinsp;4), were detected by comparing the level and sequence of miRNAs in the mutated vs. non-mutated samples. However, due to the dilution of mutants with the wild-type allele, contamination of cancer samples with normal tissues, and cancer heterogeneity, the sensitivity of the latter approach was low. Therefore, the mutations identified here likely reflect only a fraction of the functionally relevant events and should be interpreted as representative mutations rather than an exhaustive set.\u003c/p\u003e\n\u003cp\u003eIn none of the analyses did we observe that mutations affecting the miRNA genes clustered in any subregion of the gene, which could indicate the particular importance of this region for the tested functionalities. However, interesting observations may include the mutation n.84C\u0026thinsp;\u0026gt;\u0026thinsp;T in \u003cem\u003eMIR365B\u003c/em\u003e, which changes the stability of one side of miRNA duplexes and thus drastically reverses the miRNA strand balance. This observation is in agreement with the notion that the miRNA strand with the less stable 5p-end is more favorably loaded into the RISC, becoming a more stable mature miRNA\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e54\u003c/span\u003e\u003c/sup\u003e. Finally, we showed that most miRNA gene mutations decrease the stability of miRNA precursors. We also showed that mutations that we detected as having a functional effect cause, on average, a greater decrease in miRNA precursor stability. This association suggests that at least some mutation effects are expressed via changes in the miRNA precursor structure and confirms the importance of the structure and its stability for the proper functioning of miRNA precursors.\u003c/p\u003e\n\u003cp\u003eAlthough most tested mutations, including those affecting the molecular functionality of miRNA genes, are randomly occurring neutral variants, single ones may still play a role in driving cancer development. This applies, in particular, to mutations in cancer-related miRNA genes, such as those annotated in the Cancer miRNA Census\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e55\u003c/span\u003e\u003c/sup\u003e, including 28 miRNA genes annotated in CMC, in which mutations affecting miRNA gene functioning were detected. One possible example of such a gene is \u003cem\u003eMIR142\u003c/em\u003e, which is recurrently mutated in various hematological malignancies, including AML, CLL, DLBCL, FL, and other types of B-cell lymphomas, and was found as the most frequently mutated miRNA gene in any cancer\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e7\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e68\u003c/span\u003e\u003c/sup\u003e. \u003cem\u003eMIR142\u003c/em\u003e is highly expressed in lymphoid blood cells, and it has been shown that miR-142 (particularly miR-142-3p) plays an important role in hematopoiesis, regulating the development and function of different hematologic lineages as well as in hematological malignancies\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e69\u003c/span\u003e\u003c/sup\u003e. Functional studies have focused mainly on a few mutations in the miR-142-3p seed. However, these studies, conducted in various cellular and animal models, demonstrated that the mutations, through ineffective regulation of miR-142-3p targets, including \u003cem\u003eASH1L\u003c/em\u003e, increase the level of HOXA9/A10, resulting in aberrant hematopoietic differentiation\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e17\u003c/span\u003e\u003c/sup\u003e. This promotes myeloid and suppresses lymphoid lineages, ultimately leading to leukemic transformation and AML. Mutations synergize with mutations in \u003cem\u003eIDH2\u003c/em\u003e\u003csup\u003e\u003cem\u003e59\u003c/em\u003e\u003c/sup\u003e. A gain-of-function effect for one of the mutations was also suggested\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e70\u003c/span\u003e\u003c/sup\u003e. The functional effects of the seed mutations were observed despite their small impact on the level of miR-142-3p, suggesting that the mutations act predominantly by affecting the seed sequence and thus target recognition. This finding is consistent with the lack of effects of n.59T\u0026thinsp;\u0026gt;\u0026thinsp;C (two samples) and n.58G\u0026thinsp;\u0026gt;\u0026thinsp;C on the miR-142-3p level observed in our study (Experiment_1; Fig. \u003cspan class=\"InternalRef\"\u003e2\u003c/span\u003eA). A recent study published during the preparation of this manuscript also demonstrated the molecular effects of mutations located outside the 3p seed in different parts of \u003cem\u003eMIR142\u003c/em\u003e overexpressed in HEK293 cells\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e19\u003c/span\u003e\u003c/sup\u003e. To date, however, \u003cem\u003eMIR142\u003c/em\u003e mutations have been studied only in artificially generated cellular or mouse models. Here, for the first time, we demonstrated the molecular effects of \u003cem\u003eMIR142\u003c/em\u003e mutations in relevant cancer samples, predominantly LAML and DLBC. We found that two mutations, n.55A\u0026thinsp;\u0026gt;\u0026thinsp;G and n.16C\u0026thinsp;\u0026gt;\u0026thinsp;A, severely affected the levels of miR-142-3p and miR-142-5p (Experiment_1; Supplementary Table S2). Additionally, three \u003cem\u003eMIR142\u003c/em\u003e mutations significantly affected the distribution of miR-142-3p isomiRs (n.55A\u0026thinsp;\u0026gt;\u0026thinsp;G and n.59T\u0026thinsp;\u0026gt;\u0026thinsp;C in Experiment_2 and n.85G\u0026thinsp;\u0026gt;\u0026thinsp;A in Experiment_4). Surprisingly, despite being located in different parts of the precursor, all these mutations induced similar effects on the isomiR profile, leading to a decrease in +\u0026thinsp;1|n isomiRs (predominant in wild-type) in favor of 0|n isomiRs. A similar (although less profound) shift in the isomiR profile was also triggered by two other mutations in \u003cem\u003eMIR142\u003c/em\u003e (n.58G\u0026thinsp;\u0026gt;\u0026thinsp;C and n.59T\u0026thinsp;\u0026gt;\u0026thinsp;C). Thus, all of these mutations affect the miR-142-3p seed by shifting the miRNA; additionally, 3 of the mutations located in the 3p seed also directly alter its sequence.\u003c/p\u003e\n\u003cp\u003eAnother well-known cancer-related miRNA gene is \u003cem\u003eMIR205\u003c/em\u003e, which is important for cancer development. It has been assigned both an oncogenic role and a tumor suppressor role, depending on the tissue/cancer type\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e71\u003c/span\u003e,\u003cspan class=\"CitationRef\"\u003e72\u003c/span\u003e\u003c/sup\u003e. Analysis of TCGA datasets revealed enrichment of mutations in \u003cem\u003eMIR205\u003c/em\u003e\u003csup\u003e\u003cem\u003e7\u003c/em\u003e\u003c/sup\u003e in different solid cancers, particularly in melanoma; however, their impact on miRNA functioning has never been studied. Moreover, a reduced level of miR-205-5p is associated with shorter survival in melanoma patients\u003csup\u003e\u003cspan class=\"CitationRef\"\u003e73\u003c/span\u003e\u003c/sup\u003e. In this study, we found that n.35C\u0026thinsp;\u0026gt;\u0026thinsp;T, located in the 5p-arm of \u003cem\u003eMIR205\u003c/em\u003e (Experiment_1), decreases the level of miR-205 in SKCM (melanoma). Two other mutations, n.30C\u0026thinsp;\u0026gt;\u0026thinsp;T and n.94G\u0026thinsp;\u0026gt;\u0026thinsp;C, identified in BLCA and localized in 5p- and 3p-flanks, respectively, affect the isomiR profiles of miR-205-5p (Experiment_4; Supplementary Figure S5). Finally, n.97-98dupGG, located in the 3p-flank (Experiment_5, Fig. \u003cspan class=\"InternalRef\"\u003e6\u003c/span\u003eC), affects the structure of miRNA precursors, resulting in a shift in DROSHA and DICER1 cleavage sites and reversing the 5p/3p strand balance, reducing the fraction of miR-205-5p (dominant mature strand) in favor of miR-205-3p (passenger strand).\u003c/p\u003e\n\u003cp\u003eTwo more examples of mutations located in cancer-related genes are n.21-22delTA in \u003cem\u003eMIR10B\u003c/em\u003e\u003csup\u003e\u003cem\u003e\u003cspan class=\"CitationRef\"\u003e74\u003c/span\u003e\u003c/em\u003e\u003c/sup\u003e and n.-6C\u0026thinsp;\u0026gt;\u0026thinsp;G in \u003cem\u003eMIR21\u003c/em\u003e\u003csup\u003e\u003cem\u003e75\u003c/em\u003e\u003c/sup\u003e. Both mutations led to significant changes in isomiR profiles (Experiment_4; Fig. \u003cspan class=\"InternalRef\"\u003e5\u003c/span\u003e). Mutation in \u003cem\u003eMIR10B\u003c/em\u003e increases the +\u0026thinsp;1|n and +\u0026thinsp;2|n fractions of miR-10b-5p, changing its canonical seed. In contrast, the mutation in \u003cem\u003eMIR21\u003c/em\u003e affects miR-21-3p, decreasing the fraction of 0|+ in favor of 0|0 and non-templated isomiRs.\u003c/p\u003e\n\u003cp\u003eIn summary, to the best of our knowledge, we performed the broadest systematic analysis of the effects of mutations in miRNA genes, in which we identified 87 mutations that significantly affect different miRNA gene functionalities, including miRNA levels, isomiR profiles, and strand balance. For the first time, we studied the effects of mutations in natural genomic contexts, not in artificially generated models. Our results improve the understanding of the impact of genetic variants on miRNA biogenesis and may help develop tools for predicting the significance of genetic variants in miRNA genes. Generally, our results (the first approach) indicate that most miRNA gene mutations, not only those located in seeds, affect the proper functioning of miRNA genes and should therefore be considered likely to reveal deleterious variants in genetic analyses.\u003c/p\u003e"},{"header":"Declarations","content":"\u003cp\u003e\u003cstrong\u003eAuthors\u0026rsquo; contribution\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003e\u003cu\u003eMagdalena Machowska\u003c/u\u003e: Conceptualization, Methodology, Validation, Formal analysis, Resources, Data Curation, Visualization, Writing - Original Draft, Writing - Review \u0026amp; Editing. \u003cu\u003eNatalia Szostak\u003c/u\u003e: Methodology, Software, Validation, Formal analysis, Resources, Data Curation, Writing - Review \u0026amp; Editing. \u003cu\u003eAdrian Tire\u003c/u\u003e: Formal analysis, Visualization, Writing - Review \u0026amp; Editing. \u003cu\u003eWladyslaw Wegorek\u003c/u\u003e: Formal analysis, Visualization, Writing - Review \u0026amp; Editing. \u003cu\u003eMalwina Suszynska\u003c/u\u003e:\u0026nbsp;Formal analysis, Writing - Review \u0026amp; Editing. \u003cu\u003eArkadiusz Kajdasz\u003c/u\u003e:\u0026nbsp;Software, Formal analysis, Writing - Review \u0026amp; Editing. \u003cu\u003ePaulina Galka-Marciniak\u003c/u\u003e:\u0026nbsp;Formal analysis, Funding acquisition, Writing - Review \u0026amp; Editing. \u003cu\u003eAnna Philips\u003c/u\u003e:\u0026nbsp;Resources, Writing - Review \u0026amp; Editing. \u003cu\u003ePiotr Kozlowski\u003c/u\u003e:\u0026nbsp;Conceptualization, Methodology, Validation, Formal analysis, Resources, Data Curation, Visualization, Supervision, Project administration, Funding acquisition, Writing - Original Draft, Writing - Review \u0026amp; Editing.\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eAcknowledgements\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe study was supported by grants from the Polish National Science Centre (2020/39/B/NZ5/01970 and 2016/22/A/NZ2/00184 awarded to Piotr Kozlowski, and 2020/39/D/NZ2/03106 awarded to Paulina Galka-Marciniak).\u003c/p\u003e\n\u003cp\u003e\u003cstrong\u003eCompeting interests\u003c/strong\u003e\u003c/p\u003e\n\u003cp\u003eThe authors have no conflicts of interest to declare.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eGebert LFR, MacRae IJ (2019) Regulation of microRNA function in animals. Nat Rev Mol Cell Biol 20:21\u0026ndash;37\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eO\u0026rsquo;Brien J, Hayder H, Zayed Y, Peng C (2018) Overview of MicroRNA Biogenesis, Mechanisms of Actions, and Circulation. Front Endocrinol (Lausanne) 9\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSaliminejad K, Khorshid K, Soleymani Fard HR, S., Ghaffari SH (2019) An overview of microRNAs: Biology, functions, therapeutics, and analysis methods. J Cell Physiol 234:5451\u0026ndash;5465\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAli Syeda Z, Langden SSS, Munkhzul C, Lee M, Song SJ (2020) Regulatory Mechanism of MicroRNA Expression in Cancer. Int J Mol Sci 21:1723\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKozomara A, Birgaoanu M, Griffiths-Jones S (2019) miRBase: from microRNA sequences to function. Nucleic Acids Res 47:D155\u0026ndash;D162\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFromm B et al (2022) MirGeneDB 2.1: toward a complete sampling of all major animal phyla. Nucleic Acids Res 50:D204\u0026ndash;D210\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eUrbanek-Trzeciak MO et al (2020) Pan-cancer analysis of somatic mutations in miRNA genes. EBioMedicine 61:103051\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHan M, Zheng Y (2013) Comprehensive Analysis of Single Nucleotide Polymorphisms in Human MicroRNAs. PLoS ONE 8:e78028\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMachowska M, Galka-Marciniak P, Kozlowski P (2022) Consequences of genetic variants in miRNA genes. Comput Struct Biotechnol J 20:6443\u0026ndash;6457\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eUrbanek-Trzeciak MO, Kozlowski P, Galka-Marciniak P, miRMut (2022) Annotation of mutations in miRNA genes from human whole-exome or whole-genome sequencing. STAR Protoc 3:101023\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMenc\u0026iacute;a \u0026Aacute; et al (2009) Mutations in the seed region of human miR-96 are responsible for nonsyndromic progressive hearing loss. Nat Genet 41:609\u0026ndash;613\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHughes AE et al (2011) Mutation Altering the miR-184 Seed Region Causes Familial Keratoconus with Cataract. Am J Hum Genet 89:628\u0026ndash;633\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eConte I et al (2015) MiR-204 is responsible for inherited retinal dystrophy associated with ocular coloboma. \u003cem\u003eProceedings of the National Academy of Sciences\u003c/em\u003e 112\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGrigelioniene G et al (2019) Gain-of-function mutation of microRNA-140 in human skeletal dysplasia. Nat Med 25:583\u0026ndash;590\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCalin GA et al (2002) Frequent deletions and down-regulation of micro- RNA genes miR15 and miR16 at 13q14 in chronic lymphocytic leukemia. \u003cem\u003eProceedings of the National Academy of Sciences\u003c/em\u003e 99, 15524\u0026ndash;15529\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCimmino A et al (2005) miR-15 and miR-16 induce apoptosis by targeting BCL2. \u003cem\u003eProceedings of the National Academy of Sciences\u003c/em\u003e 102, 13944\u0026ndash;13949\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTrissal MC et al (2018) MIR142 Loss-of-Function Mutations Derepress ASH1L to Increase HOXA Gene Expression and Promote Leukemogenesis. Cancer Res 78:3510\u0026ndash;3521\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKwanhian W et al (2012) Micro RNA -142 is mutated in about 20% of diffuse large B‐cell lymphoma. Cancer Med 1:141\u0026ndash;155\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBayraktar R et al (2025) The mutational landscape and functional effects of noncoding ultraconserved elements in human cancers. Sci Adv 11\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSun G et al (2009) SNPs in human miRNA genes affect biogenesis and function. RNA 15:1640\u0026ndash;1651\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLeit\u0026atilde;o AL, Enguita FJ (2022) A Structural View of miRNA Biogenesis and Function. Noncoding RNA 8:10\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTreiber T, Treiber N, Meister G (2019) Regulation of microRNA biogenesis and its crosstalk with other cellular pathways. Nat Rev Mol Cell Biol 20:5\u0026ndash;20\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShang R, Lee S, Senavirathne G, Lai EC (2023) microRNAs in action: biogenesis, function and regulation. Nat Rev Genet 24:816\u0026ndash;833\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBaek SC et al (2024) Structural atlas of human primary microRNAs generated by SHAPE-MaP. Mol Cell 84:1158\u0026ndash;1172e6\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRoden C et al (2017) Novel determinants of mammalian primary microRNA processing revealed by systematic evaluation of hairpin-containing transcripts and human genetic variation. Genome Res 27:374\u0026ndash;384\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHalvorsen M, Martin JS, Broadaway S, Laederach A (2010) Disease-associated mutations that alter the RNA structural ensemble. PLoS Genet 6:e1001074\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOmariba G et al (2020) Genome-Wide Analysis of MicroRNA-related Single Nucleotide Polymorphisms (SNPs) in Mouse Genome. Sci Rep 10:5789\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGong J et al (2012) Genome-wide identification of SNPs in microRNA genes and the SNP effects on microRNA target binding and biogenesis. Hum Mutat 33:254\u0026ndash;263\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFang W, Bartel DP (2015) The Menu of Features that Define Primary MicroRNAs and Enable De Novo Design of MicroRNA Genes. Mol Cell 60:131\u0026ndash;145\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRice GM, Shivashankar V, Ma EJ, Baryza JL, Nutiu R (2020) Functional Atlas of Primary miRNA Maturation by the Microprocessor. Mol Cell 80:892\u0026ndash;902e4\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAuyeung VC, Ulitsky I, McGeary SE, Bartel DP (2013) Beyond secondary structure: primary-sequence determinants license pri-miRNA hairpins for processing. Cell 152:844\u0026ndash;858\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGalka-Marciniak P et al (2016) siRNA release from pri-miRNA scaffolds is controlled by the sequence and structure of RNA. Biochim Biophys Acta 1859:639\u0026ndash;649\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStarega-Roslan J, Witkos T, Galka-Marciniak P, Krzyzosiak W (2015) Sequence Features of Drosha and Dicer Cleavage Sites Affect the Complexity of IsomiRs. Int J Mol Sci 16:8110\u0026ndash;8127\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNguyen TD, Trinh TA, Bao S, Nguyen TA (2022) Secondary structure RNA elements control the cleavage activity of DICER. Nat Commun 13:2138\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKim K et al (2021) A quantitative map of human primary microRNA processing sites. Mol Cell 81:3422\u0026ndash;3439e11\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLoher P et al (2021) IsoMiRmap: fast, deterministic and exhaustive mining of isomiRs from short RNA-seq datasets. Bioinformatics 37:1828\u0026ndash;1838\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHoadley KA et al (2018) Cell-of-Origin Patterns Dominate the Molecular Classification of 10,000 Tumors from 33 Types of Cancer. Cell 173:291\u0026ndash;304e6\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eIbing S et al (2021) On the impact of batch effect correction in TCGA isomiR expression data. NAR Cancer 3\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eColaprico A et al (2016) TCGAbiolinks: an R/Bioconductor package for integrative analysis of TCGA data. Nucleic Acids Res 44:e71\u0026ndash;e71\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLeek JT, Johnson WE, Parker HS, Jaffe AE, Storey JD (2012) The sva package for removing batch effects and other unwanted variation in high-throughput experiments. Bioinformatics 28:882\u0026ndash;883\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRitchie ME et al (2015) limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res 43:e47\u0026ndash;e47\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZuker M (2003) Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res 31:3406\u0026ndash;3415\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDarty K, Denise A, Ponty YVARNA (2009) Interactive drawing and editing of the RNA secondary structure. Bioinformatics 25:1974\u0026ndash;1975\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGalka-Marciniak P et al (2019) Somatic Mutations in miRNA Genes in Lung Cancer\u0026mdash;Potential Functional Consequences of Non-Coding Sequence Variants. Cancers (Basel) 11:793\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePopenda M et al (2012) Automated 3D structure composition for large RNAs. Nucleic Acids Res 40:e112\u0026ndash;e112\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eThe PyMOL Molecular Graphics System, Version 3.0 Schr\u0026ouml;dinger, LLC.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLoher P, Londin ER, Rigoutsos I (2014) IsomiR expression profiles in human lymphoblastoid cell lines exhibit population and gender dependencies. Oncotarget 5:8790\u0026ndash;8802\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRea LM, Parker RA (1992) Designing and Conducting Survey Research. Jossey-Bass, San Francisco\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLove MI, Huber W, Anders S (2014) Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol 15:550\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVaret H, Brillet-Gu\u0026eacute;guen L, Copp\u0026eacute;e J-Y, Dillies M-A, SARTools (2016) A DESeq2- and EdgeR-Based R Pipeline for Comprehensive Differential Analysis of RNA-Seq Data. PLoS ONE 11:e0157022\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAgarwal V, Bell GW, Nam J-W (2015) \u0026amp; Bartel, D. P. Predicting effective microRNA target sites in mammalian mRNAs. \u003cem\u003eElife\u003c/em\u003e 4\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChen Y, Wang X (2020) miRDB: an online database for prediction of functional microRNA targets. Nucleic Acids Res 48:D127\u0026ndash;D131\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLee Y-Y, Kim H, Kim V (2023) N. Sequence determinant of small RNA production by DICER. Nature 615:323\u0026ndash;330\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMedley JC, Panzade G, Zinovyeva A (2021) Y. microRNA strand selection: Unwinding the rules. WIREs RNA 12\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSuszynska M et al (2024) Cancer miRNA Census \u0026ndash; a list of cancer-related miRNA genes. Nucleic Acids Res 52:1628\u0026ndash;1644\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDietlein F et al (2022) Genome-wide analysis of somatic noncoding mutation patterns in cancer. Sci (1979) 376\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHornsh\u0026oslash;j H et al (2018) Pan-cancer screen for mutations in non-coding elements with conservation and cancer specificity reveals correlations with expression and survival. NPJ Genom Med 3:1\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eRheinbay E et al (2020) Analyses of non-coding somatic drivers in 2,658 cancer whole genomes. Nature 578:102\u0026ndash;111\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMarshall A et al (2020) Mir142 loss unlocks IDH2R140-dependent leukemogenesis through antagonistic regulation of HOX genes. Sci Rep 10:19390\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSold\u0026agrave; G et al (2012) A novel mutation within the MIR96 gene causes non-syndromic inherited hearing loss in an Italian family by altering pre-miRNA processing. Hum Mol Genet 21:577\u0026ndash;585\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDuan R, Pak C, Jin P (2007) Single nucleotide polymorphism associated with mature miR-125a alters the processing of pri-miRNA. Hum Mol Genet 16:1124\u0026ndash;1131\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDing W et al (2018) A polymorphism rs3746444 within the pre-miR‐499 alters the maturation of miR‐499‐5p and its antiapoptotic function. J Cell Mol Med 22:5418\u0026ndash;5428\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOhanian M, Humphreys DT, Anderson E, Preiss T, Fatkin D (2013) A heterozygous variant in the human cardiac miR-133 gene, MIR133A2, alters miRNA duplex processing and strand abundance. BMC Genet 14:18\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhao D et al (2022) Exosomal miR-1304-3p promotes breast cancer progression in African Americans by activating cancer-associated adipocytes. Nat Commun 13:7734\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLechner J et al (2013) Mutational Analysis of MIR184 in Sporadic Keratoconus and Myopia. Invest Opthalmology Visual Sci 54:5266\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eCalin GA et al (2005) A MicroRNA Signature Associated with Prognosis and Progression in Chronic Lymphocytic Leukemia. N Engl J Med 353:1793\u0026ndash;1801\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eIwai N, Naraba H (2005) Polymorphisms in human pre-miRNAs. Biochem Biophys Res Commun 331:1439\u0026ndash;1444\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGalka-Marciniak P et al (2022) Mutations in the miR-142 gene are not common in myeloproliferative neoplasms. Sci Rep 12:10924\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHuang W, Paul D, Calin GA, Bayraktar R (2023) miR-142: A Master Regulator in Hematological Malignancies and Therapeutic Opportunities. Cells 13:84\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKawano S et al (2023) A gain-of‐function mutation in microRNA 142 is sufficient to cause the development of T‐cell leukemia in mice. Cancer Sci 114:2821\u0026ndash;2834\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChauhan N, Dhasmana A, Jaggi M, Chauhan SC, Yallapu MM (2020) miR-205: A Potential Biomedicine for Cancer Therapy. Cells 9:1957\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFerrari E, Gandellini P (2020) Unveiling the ups and downs of miR-205 in physiology and cancer: transcriptional and post-transcriptional mechanisms. Cell Death Dis 11:980\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eS\u0026aacute;nchez-Sendra B et al (2018) Downregulation of intratumoral expression of miR-205, miR-200c and miR-125b in primary human cutaneous melanomas predicts shorter survival. Sci Rep 8:17076\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSheedy P, Medarova Z (2018) The fundamental role of miR-10b in metastatic cancer. Am J Cancer Res 8:1674\u0026ndash;1688\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBautista-S\u0026aacute;nchez D et al (2020) The Promising Role of miR-21 as a Cancer Biomarker and Its Importance in RNA-Based Therapeutics. Mol Ther Nucleic Acids 20:409\u0026ndash;420\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":true,"hideJournal":true,"highlight":"","institution":"","isAcceptedByJournal":false,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true},"keywords":"","lastPublishedDoi":"10.21203/rs.3.rs-7029847/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-7029847/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003cp\u003eA growing number of mutations are being identified in the noncoding genome, including miRNA genes. However, little is known about the consequences of these mutations and how harmful they are to the functioning of miRNA genes. To evaluate the effects of miRNA gene mutations, we took advantage of a large collection of somatic mutations identified in miRNA genes in \u0026gt;\u0026thinsp;10K TCGA cancer samples and compared them with the corresponding miRNA-seq data. Using different analytical approaches and highly rigorous statistical criteria, we identified many mutations (n\u0026thinsp;=\u0026thinsp;87) that affect the level of mature miRNAs (predominantly decreasing), isomiR profiles (precision of DROSHA/DICER1 cleavage), and/or 5p/3p miRNA strand balance. Taken together, the analysis revealed that most miRNA gene mutations, not only those in the seed, may be deleterious for the proper functioning of miRNA genes. We also showed that most miRNA gene mutations destabilize the structure of miRNA precursors and that mutations identified as deleterious are associated with a stronger destabilizing effect. Moreover, although most cancer somatic mutations are randomly occurring neutral variants, some mutations that alter the function of well-known cancer-related miRNA genes, such as \u003cem\u003eMIR21\u003c/em\u003e, \u003cem\u003eMIR142\u003c/em\u003e, or \u003cem\u003eMIR205\u003c/em\u003e, might be functional variants in cancer.\u003c/p\u003e","manuscriptTitle":"miRNA gene mutations commonly disrupt the proper functioning of miRNA genes","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-07-21 08:54:19","doi":"10.21203/rs.3.rs-7029847/v1","editorialEvents":[{"type":"communityComments","content":0}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"researchsquare","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":true,"externalIdentity":"","sideBox":"","snPcode":"","submissionUrl":"/submission","title":"Research Square","twitterHandle":"researchsquare","acdcEnabled":true,"dfaEnabled":false,"editorialSystem":"","reportingPortfolio":"","inReviewEnabled":false,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"710f6b75-b66f-45e7-bf80-1377731ce0be","owner":[],"postedDate":"July 21st, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"posted","subjectAreas":[{"id":51488279,"name":"Biological sciences/Genetics/Functional genomics"},{"id":51488280,"name":"Biological sciences/Molecular biology/Non-coding RNAs/miRNAs"},{"id":51488281,"name":"Biological sciences/Genetics/Gene regulation"},{"id":51488282,"name":"Biological sciences/Genetics/Mutation"},{"id":51488283,"name":"Biological sciences/Computational biology and bioinformatics"}],"tags":[],"updatedAt":"2025-08-04T20:25:54+00:00","versionOfRecord":[],"versionCreatedAt":"2025-07-21 08:54:19","video":"","vorDoi":"","vorDoiUrl":"","workflowStages":[]},"version":"v1","identity":"rs-7029847","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-7029847","identity":"rs-7029847","version":["v1"]},"buildId":"8U1c8b4HqxoKbykW_rLl7","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}