{"id":628000,"date":"2023-04-11T11:56:45","date_gmt":"2023-04-11T16:56:45","guid":{"rendered":"https:\/\/news.sellorbuyhomefast.com\/index.php\/2023\/04\/11\/rip-pen-seq-identifies-a-class-of-kink-turn-rnas-as-splicing-regulators\/"},"modified":"2023-04-11T11:56:45","modified_gmt":"2023-04-11T16:56:45","slug":"rip-pen-seq-identifies-a-class-of-kink-turn-rnas-as-splicing-regulators","status":"publish","type":"post","link":"https:\/\/newsycanuse.com\/index.php\/2023\/04\/11\/rip-pen-seq-identifies-a-class-of-kink-turn-rnas-as-splicing-regulators\/","title":{"rendered":"RIP-PEN-seq identifies a class of kink-turn RNAs as splicing regulators"},"content":{"rendered":"

Science & Nature <\/p>\n

\n
\n

Main<\/h2>\n
\n

RNA structural motifs confer diversity of RNA architecture to regulate various biological processes1<\/a><\/sup>. The kink-turn (K-turn) is the most prevalent three-dimensional (3D) RNA structural motif in messenger RNAs (mRNAs) and non-coding RNAs (ncRNAs)2<\/a>,3<\/a>,4<\/a><\/sup>. The K-turn structures are typically characterized by a canonical stem (C-stem) followed by an asymmetric 3-nucleotide (3-nt) bulge and G\u2022A and A\u2022G base pairs on the 3\u2032 non-canonical stem (NC-stem)2<\/a>,3<\/a><\/sup>. The naturally occurring K-turn often serves as a specific binding site for RNA-binding proteins (RBPs), including the protein L7Ae5<\/a><\/sup> in archaea and the other homologues YbxF, YlxQ6<\/a><\/sup> and 15.5-kDa protein (15.5K)7<\/a><\/sup>. Notably, the unique conformational property of the K-turn, as well as its binding proteins, enables the K-turn to play pivotal roles in maintaining the proper structure and biological function of various types of RNAs2<\/a>,3<\/a>,4<\/a><\/sup>.<\/p>\n

The K-turn structural motifs are recurrently observed in various types of ncRNA families and play important roles in RNA metabolism2<\/a>,3<\/a>,4<\/a><\/sup>. In eukaryotes, the well-known K-turn ncRNAs are the box C\/D small nucleolar RNAs (snoRNAs) and U4\/U4atac small nuclear RNAs (snRNAs). The C\/D boxes are located at the box C\/D RNA terminal pair and form the terminal K-turn motif8<\/a><\/sup>. Recognition of the K-turn structural motif of C\/D box RNAs by the 15.5K protein initiates small nucleolar ribonucleoprotein (snoRNP) assembly to carry out site-specific 2\u2032-O<\/i>-methylation of ribosomal RNAs (rRNAs) and snRNAs9<\/a><\/sup>, which is essential for the functional fidelity of ribosomes and gene expression10<\/a>,11<\/a><\/sup>. The 15.5K protein also binds to the identical K-turn structural motifs in the 5\u2032 stem-loop of U4 and U4atac, promoting the assembly of the major and minor spliceosomes12<\/a>,13<\/a><\/sup>. Mutations in the U4atac K-turn structural motif have been shown to impair the binding of the 15.5K protein14<\/a><\/sup> and cause various diseases, such as the Taybi\u2013Linder syndrome (TALS\/MOPD1)15<\/a>,16<\/a><\/sup>, Roifman syndrome (RFMN)17<\/a><\/sup> and Lowry\u2013Wood syndrome (LWS)18<\/a><\/sup>. However, the prevalence, mechanism and function of K-turn structures in the transcriptome remain largely unknown.<\/p>\n

In this study, we developed a method for RNA immunoprecipitation coupled with sequencing of paired ends of ncRNAs (RIP-PEN-seq) to identify the full-length sequences of ncRNAs bound by the 15.5K protein and discover a class of backward K-turn RNAs (bktRNAs; Fig. 1a<\/a>) with two consensus motifs. In contrast to all well-known K-turn RNAs (ktRNAs) in humans, whose internal loops are located at the 5\u2032 ends of the RNAs (forward ktRNAs (fktRNAs); Fig. 1a<\/a>), the internal loops of the ktRNAs identified in this study are located at the 3\u2032 termini; we, therefore, named them bktRNAs (Fig. 1a<\/a>). These bktRNAs have a strong preference for a 5\u2032 motif (CUGA) and a 3\u2032 motif (UGAUG) at the fixed terminal position of RNA. Moreover, we identified a functional bktRNA1 that interacts with the 15.5K and fibrillarin (FBL) proteins to regulate RNA modification of U12 snRNA. Loss of bktRNA1 results in aberrant RNA splicing of hundreds of U12-type introns. In addition, we reveal bktRNA1 and that 80% of randomly selected bktRNAs are involved in the local regulation of intron splicing. Together, our results reveal that bktRNAs form a previously undescribed class of post-transcriptional regulators.<\/p>\n

\n
Fig. 1: A class of ncRNAs with a backward K-turn structure.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n

a<\/b>, Secondary structure of the consensus fktRNA and bktRNA. X represents any nucleotide. The helix 5\u2032 to the loop is called the C (canonical) stem (C-stem), and the G\u2022A and A\u2022G base pairs at the 3\u2032 side are called the NC-stem. There is a 3-nt internal loop between the C-stem and the NC-stem. In the upper panel, the internal loop and the C-stem located at the 5\u2032 end of the RNA compose the fktRNA, whereas the internal loop and the C-stem located at the 3\u2032 terminus compose the bktRNA. b<\/b>, Procedure for the construction of RIP-PEN-seq libraries. c<\/b>,d<\/b>, Histogram showing the numbers of bktRNAs with different distances between the backward K-turn structure and the bktRNA 5\u2032 end (c<\/b>) or 3\u2032 end (d<\/b>). e<\/b>, Enriched motifs identified in bktRNAs by MEME software48<\/a><\/sup>. Two significant motifs were identified in the bktRNAs. One motif (CUGA) is near the 5\u2032 end of the bktRNA, and one motif (UGAUG) is near the 3\u2032 end of the bktRNA. The motif P<\/i> value is estimated based on the E<\/i> value of the log-likelihood ratio. f<\/b>, Gene model of bktRNA identified from RIP-PEN-seq data. g<\/b>, Violin plots displaying the SHAPE reactivity across the backward K-turn structure (including 5\u2032 NC stem, 5\u2032 C stem, Loop, 3\u2032 C stem, internal (internal loop) and 3\u2032 NC stem), averaged across all bktRNAs (n<\/i>\u2009=\u200922) in the RIP-PEN-SHAPE-MaP data. The box plots indicate the median and the upper and lower quartiles. h<\/b>,i<\/b>, The predicted secondary structure of the backward K-turn (left panel) and SHAPE reactivity signal (right panel) on bktRNA1 (h<\/b>) and bktRNA2 (i<\/b>). The 5\u2032 and 3\u2032 motifs are indicated with black boxes in the structure figures. The NC-stem and C-stem are marked with black and red underlines in the bar plot, respectively. The SHAPE reactivity is calculated from merged n<\/i>\u2009=\u20094 biological replicates. RT, reverse transcription.<\/p>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n<\/div>\n<\/div>\n

\n

Results<\/h2>\n
\n

bktRNAs with consensus motifs identified by RIP-PEN-seq<\/h3>\n

To identify ncRNAs (20\u2013500\u2009nt) with K-turn structural motifs, we developed an efficient RNA cloning scheme called RIP-PEN-seq that uses dual RNA adapters and size selection, and a series of advanced experimental strategies to capture both ends of any ncRNA bound by the 15.5K protein (Fig. 1b<\/a>, Extended Data Fig. 1a<\/a> and Methods<\/a>). We applied RIP-PEN-seq to HEK293T cells expressing FLAG-tagged 15.5K (Extended Data Fig. 1b,c<\/a>), and a total of approximately 520 million paired-end reads were obtained, of which approximately 282 million could be uniquely mapped to the human genome. Notably, more than 45% of the uniquely aligned reads mapped to previously defined fktRNAs. The start and end sites of these known fktRNAs could be precisely identified at single-nucleotide resolution (Extended Data Fig. 1d,e<\/a>). For example, the precise 5\u2032 and 3\u2032 ends of 10 known C\/D box fktRNAs with K-turn motifs located within the introns of GAS5<\/i> genes were identified by our RIP-PEN-seq method (Extended Data Fig. 1f<\/a>). Together, these results suggest that our RIP-PEN-seq approach not only shows high specificity and accuracy in enriching ktRNAs but can also capture their full-length sequences.<\/p>\n

Although almost all known ktRNAs are fktRNAs in humans, as the 15.5K protein has no preference for binding forward or backward K-turn structural motifs2<\/a>,19<\/a><\/sup>, we hypothesized that the 15.5K RIP-PEN-seq data may contain hidden bktRNAs. We then developed a computational tool, called kturnSeeker, to identify potential ktRNAs based on sequence and secondary structure information from the RIP-PEN-seq data (Extended Data Fig. 1g,h<\/a>). As a result, we identified 605 previously undescribed fktRNAs as well as 118 canonical fktRNAs (Extended Data Fig. 1i<\/a> and Supplementary Table 1<\/a>). Notably, most previously undescribed fktRNAs are located within introns (Extended Data Fig. 1i<\/a> and Supplementary Table 1<\/a>). Furthermore, after removing fktRNAs and known annotations and keeping the ktRNAs present in at least two sequencing libraries, we discovered 84 ktRNAs (alignment to 93 genomic loci) with backward K-turn structural motifs, which we named bktRNA1 to bktRNA84 (Supplementary Table 1<\/a>).<\/p>\n

Intriguingly, the K-turn structural motifs of most bktRNAs are located at 4\u2009nt and 2\u2009nt from the RNA 5\u2032 end and 3\u2032 end, respectively (Fig. 1c,d<\/a>). Moreover, we observed highly substantial enrichment of a CUGA motif often 4\u2009nt downstream of the 5\u2032 end and a UGAUG motif 2\u2009nt upstream of the 3\u2032 end by de novo motif enrichment analysis (Fig. 1e<\/a> and Extended Data Fig. 2a<\/a>). Notably, all these enrichment consensus motifs were located within the K-turn structural motifs of bktRNAs (Fig. 1f\u2013i<\/a> and Extended Data Fig. 2b<\/a>).<\/p>\n

To further confirm that these two conserved motifs can truly adopt the backward K-turn structure in vivo, we developed 15.5K RIP-PEN-SHAPE-MaP, an approach combining SHAPE-MaP20<\/a>,21<\/a><\/sup> and our RIP-PEN-seq experiments to probe intact RNA structures of bktRNAs (Extended Data Fig. 2c<\/a> and Methods<\/a>). We first constructed reactivity scores for the known fktRNAs and found that paired nucleotides in the NC-stem and C-stems exhibited lower SHAPE reactivity and that unpaired nucleotides in the loop region (Loop) had relatively high SHAPE reactivity (Extended Data Fig. 2d,e<\/a> and Supplementary Table 2<\/a>), demonstrating the accuracy of the 15.5K RIP-PEN-SHAPE-MaP method. For example, the SHAPE reactivity signal agreed well with the structure of fktRNA SNORD102 (Extended Data Fig. 2e<\/a>). We then mapped the SHAPE reactivity signals to the predicted secondary structures of bktRNAs and found that the paired nucleotides in backward K-turn structures had lower reactivity signals than the unpaired nucleotides (Fig. 1g\u2013i<\/a> and Supplementary Table 2<\/a>). Therefore, these data confirmed that the 5\u2032 and 3\u2032 motifs identified in bktRNAs are likely to adopt the backward K-turn structure in vivo.<\/p>\n

We next applied RIP-PEN-seq to mouse Hepa1-6 cells stably expressing FLAG-tagged 15.5K (Supplementary Fig. 1a,b<\/a>), and it also detected mouse fktRNAs with high specificity and accuracy (Supplementary Fig. 1c\u2013e<\/a>). We identified 43 bktRNAs in mice that shared common genomic and motif characteristics with human bktRNAs (Supplementary Fig. 1f\u2013k<\/a> and Supplementary Table 3<\/a>). Overall, by developing the RIP-PEN-seq approach and the kturnSeeker algorithm, we uncovered a class of ktRNAs with a backward K-turn structural motif composed of two consensus sequence elements in fixed positions of RNAs in humans and mice.<\/p>\n

bktRNAs have complex folding, expression and evolution patterns<\/h3>\n

To identify bktRNAs across human tissues and cells, we further applied our kturnSeeker tool to 28 small RNA sequencing (sRNA-seq) datasets generated by our PEN-seq method (Methods<\/a>) and 234 public sRNA-seq datasets produced by the ENCODE consortium. This method revealed 379 high-confidence candidate bktRNAs (Fig. 2a<\/a>), 72 of which overlapped with the abovementioned bktRNAs bound by the 15.5K protein (~86%, 72\/84; Fig. 2a<\/a>). By integrating the RIP-PEN-seq dataset and the sRNA-seq datasets, we identified 391 bktRNAs (407 genomic loci; Fig. 2b<\/a> and Supplementary Table 4<\/a>). As described for the 15.5K-bound bktRNAs, these bktRNAs had substantial enrichment of CUGA motifs 4\u2009nt downstream of the 5\u2032 end and UGAUG motifs 2\u2009nt upstream of the 3\u2032 end (Fig. 2c,d<\/a>).<\/p>\n

\n
Fig. 2: Complex expression patterns of bktRNAs in tissues, cell lines and subcellular compartments and their evolution patterns.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n

a<\/b>, Overlap analysis of bktRNAs identified from 15.5K RIP-PEN-seq and PEN-seq experiments or published data. b<\/b>, Circos plot illustrating the genomic positions, tissue expression repertoires, expression levels and evolutionary conservation of bktRNAs. The plot legend is shown in the lower panel. Red, log2<\/sub> reads per million reads (RPM)\u2009>\u20090; blue, log2<\/sub> RPM\u2009<\u20090. c<\/b>,d<\/b>, Histogram showing the numbers of bktRNAs with different distances between the backward K-turn structure and the bktRNA 5\u2032 end (c<\/b>) or 3\u2032 end (d<\/b>). e<\/b>, Length distribution of all bktRNAs. The dashed line shows the median length of bktRNAs. f<\/b>, Numbers of bktRNAs in different genomic annotation categories. g<\/b>, Heat map showing the expression profiles of bktRNAs in various human tissues using sRNA-seq data from ENCODE. The expression levels of bktRNAs in cells were categorized into corresponding tissues. Cells\/tissues with at least 10 datasets were reserved. h<\/b>, Heat map of bktRNA expression in PEN-seq libraries. i<\/b>, Expression profiles of bktRNAs in different subcellular regions (cytoplasm, nucleolus and nucleoplasm) in HCT116 and HEK293T cells. j<\/b>, Violin plots displaying the expression intensity distribution of bktRNAs (n<\/i>\u2009=\u2009355) in different subcellular regions in HCT116 and HEK293T cells. The box plots indicate the median and the upper and lower quartiles. Cyto, cytoplasm; No, nucleolus; Np, nucleoplasm. The P<\/i> value between the two categories were determined by the two-sided Mann\u2013Whitney\u2013Wilcoxon test. k<\/b>, Simplified phylogenetic trees of human bktRNAs. Internal branches and roots, numbers of 1-1 orthologous bktRNA families for the indicated species. Tree tips, bktRNA numbers for each species.<\/p>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n

To determine the folding characteristics22<\/a>,23<\/a><\/sup> of these K-turn sequences, we first devised a nomenclature for the nucleotide positions in the backward K-turn structure according to the nomenclature for the forward K-turn structure23<\/a>,24<\/a><\/sup> (Extended Data Fig. 3a<\/a>). We further examined the distribution of 3b:3n and \u22121b:\u22121n sequences in bktRNAs (Extended Data Fig. 3b<\/a>). We found that almost all bktRNAs have a strong propensity for the 3b:3n\u2009=\u2009U:U (99.5%, 405\/407), which leads to an N3 conformation and the weak folding in response to metal ions22<\/a>,23<\/a>,24<\/a>,25<\/a><\/sup>. We also discovered that approximately 69% (281\/407) of bktRNAs have \u22121b:\u22121n\u2009=\u2009C:G or A:U, which leads to good folding in response to metal ions22<\/a>,23<\/a>,24<\/a>,25<\/a><\/sup> (Extended Data Fig. 3c<\/a>). Moreover, bktRNAs in mice showed similar folding properties (Extended Data Fig. 3d,e<\/a>). Notably, the 3b:3n\u2009=\u2009U:U identified in bktRNAs is considerably different from the distribution of ribosomal Kt-7, SAM-I riboswitch and U4 snRNAs22<\/a>,23<\/a>,24<\/a>,25<\/a><\/sup>.<\/p>\n

We next explored the m6<\/sup>A modification of bktRNAs26<\/a>,27<\/a><\/sup> by intersecting our bktRNAs with the annotated m6<\/sup>A sites from RMBase28<\/a><\/sup> and found that 12 and 1 bktRNAs contained m6<\/sup>A modifications in humans and mice (Extended Data Fig. 3f<\/a>), respectively. Moreover, neither human nor mouse bktRNAs contain m6<\/sup>A modifications at the 1n position (Extended Data Fig. 3g<\/a>). Together, these findings suggest that the sequence composition of the backward K-turn structure as well as the interaction with RBPs might contribute to the folding of bktRNAs.<\/p>\n

We further characterized bktRNAs and found that the sizes of these bktRNAs ranged from 20\u2009nt to ~400\u2009nt, with a median length of 40\u2009nt (Fig. 2e<\/a>). Notably, more than 84% of bktRNAs (342\/407) originated from the intron locus (Fig. 2f<\/a>). Intriguingly, the host genes of the bktRNAs were significantly enriched in various Gene Ontology terms, such as the actin filament-based process and regulation of cell adhesion terms (Extended Data Fig. 3h<\/a>).<\/p>\n

To explore the expression patterns of bktRNAs, we analyzed bktRNA expression levels across cells and tissues as well as subcellular localization. Numerous bktRNAs appeared to be specifically expressed across human tissues or cells (Fig. 2g,h<\/a> and Extended Data Fig. 3i<\/a>). For example, the bktRNA198<\/i> and bktRNA126<\/i> genes were substantially expressed in lung tissues (Extended Data Fig. 3i<\/a>). Moreover, we applied our PEN-seq method in subcellular RNA fractions (sub-PEN-seq) and discovered that most bktRNAs exhibited nuclear-specific localization patterns (Fig. 2i,j<\/a> and Supplementary Table 4<\/a>).<\/p>\n

We next reconstructed homologous families based on sequence similarity, and we inferred a stringent minimum evolutionary age of bktRNAs. We discovered that most (145\/391, 37.1%) bktRNA families were primate specific, but 16 (4.1%) families likely originated more than 90 million years (Myr) ago, and one (0.26%) family likely originated more than 450\u2009Myr ago (Fig. 2k<\/a>). Interestingly, we found that bktRNA1<\/i>, which is located in an intron of CWF19L1<\/i>, was the most conserved among all bktRNAs (Fig. 2k<\/a>) and is an endogenous chimeric ncRNA composed of a backward K-turn domain and an H\/ACA domain (SNORA12) (Extended Data Fig. 4a\u2013c<\/a>, Supplementary Figs. 2<\/a> and 3<\/a> and Supplementary Note 1<\/a>). Notably, the SHAPE reactivity signals from RIP-PEN-SHAPE-MaP experiments also confirmed this conservative secondary structure of bktRNA1 (Extended Data Fig. 4b,c<\/a>). Moreover, we discovered a K-turn-like structure that formed a helix\u2013internal-loop\u2013helix and did not have G\u2022A and A\u2022G base pairs in bktRNA1 (Extended Data Fig. 4b,c<\/a> and Supplementary Note 1<\/a>). In summary, we have confidently identified a large number of bktRNAs with complex folding, expression patterns and sequence conservation.<\/p>\n

U12 snRNA is the direct target of bktRNA1<\/h3>\n

To identify direct targets of bktRNAs, we performed ultraviolet cross-linking, ligation and sequencing of hybrids (CLASH) experiments for 15.5K (Fig. 3a<\/a>). By analyzing the chimeric reads in the 15.5K CLASH sequencing data (Supplementary Fig. 4a<\/a>), we found that bktRNA1 interacted with U12 snRNA and formed a duplex containing 14 perfect base pairs (Fig. 3b,c<\/a> and Supplementary Fig. 4a,b<\/a>). Moreover, PARIS data in HEK293T cells29<\/a><\/sup> also showed that bktRNA1 and U12 formed an RNA duplex in vivo (Fig. 3d<\/a> and Supplementary Fig. 5a,b<\/a>). We also found that the evolutionary conservation of the pairing region between bktRNA1 and U12 is very highly conserved across 100 vertebrates (Fig. 3e,f<\/a> and Extended Data Fig. 5a<\/a>). In addition, the fluorescence in situ hybridization (FISH) and immunofluorescence (IF) experiments confirmed that bktRNA1 and U12 snRNA showed co-localization with 15.5K proteins (Extended Data Fig. 5b,c<\/a>). These results reveal that U12 snRNA is the target of bktRNA1 bound by the 15.5K protein.<\/p>\n

\n
Fig. 3: bktRNA1 guides the FBL protein to introduce 2\u2032-O<\/i>-methylation in U12 snRNA.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n
\n

a<\/b>, Schematic of 15.5K infrared-CLASH technology. b<\/b>, 15.5K-CLASH technology was used to identify the interactions between bktRNA1 and U12 snRNA. The corresponding chimeric reads are marked with different colors; the blue parts of the chimeric reads are located within bktRNA1, and the red parts of chimeric reads are located within U12 snRNA. The curved lines connecting the two RNAs indicate the base pairings of the corresponding nucleotides. c<\/b>, Base pairing between bktRNA1 and U12 snRNA identified from CLASH data. A8 is the 2\u2032-O<\/i>-methylated site (red) in U12. d<\/b>, High-throughput sequencing data from PARIS experiments were used to identify the interactions between bktRNA1 and U12 snRNA. e<\/b>,f<\/b>, Evolutionary conservation (from human to zebrafish) of the interaction region in U12 (e<\/b>) and bktRNA1 (f<\/b>). The base-pairing regions are indicated with red dashed boxes. g<\/b>, The FBL-CLASH approach was used to identify the interactions between bktRNA1 and U12 snRNA. h<\/b>, PCR analyses of genomic DNA obtained from bktRNA1<\/i> KO cells to validate bktRNA1<\/i> biallelic deletion clones. The DNA size marker is shown in base pairs. i<\/b>, qPCR for the relative expression of bktRNA1 and U12 in bktRNA1\u2212\/\u2212<\/sup> and WT HCT116 cells. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. j<\/b>, irNorthern blotting of bktRNA1 in bktRNA1\u2212\/\u2212<\/sup> and WT HCT116 cells. U6 snRNA served as a loading control. k<\/b>, irPE analysis of U12 in bktRNA1\u2212\/\u2212<\/sup> and WT HCT116 cells. The black triangle represents the presence of decreasing amounts (0.5\u2009mM, 0.05\u2009mM and 0.01\u2009mM) of dNTPs. The representative RNA sequence is shown in the left panel, and the red nucleotides denote 2\u2032-O<\/i>-methylation sites detected by irPE. The red triangle in the right panel indicates the primer extension stop signal. l<\/b>, irNorthern blotting of bktRNA1 in bktRNA1\u2212\/\u2212<\/sup> and bktRNA1-re-expressing cells. U6 snRNA served as a loading control. m<\/b>, irPE analysis of U12 in bktRNA1\u2212\/\u2212<\/sup> and bktRNA1-re-expressing cells. Cons, conservation; RT, reverse transcription; verts, vertebrates.<\/p>\n

Source data<\/a><\/p>\n<\/div>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n

bktRNA1 is indispensable for 2\u2032-O<\/i>-methylation of U12 snRNA<\/h3>\n

As a known 2\u2032-O<\/i>-methylation site30<\/a><\/sup> is located within the complementary region between U12 snRNA and bktRNA1, we hypothesized that bktRNA1 guides the 2\u2032-O<\/i>-methyltransferase FBL to introduce 2\u2032-O<\/i>-methylation of U12 snRNA. To test this hypothesis, we first developed an infrared primer extension (irPE) method (Methods<\/a>) and confirmed that the A8 residue in U12 snRNA is 2\u2032-O<\/i>-methylated in human cells (Supplementary Fig. 6a\u2013c<\/a>). The known 2\u2032-O<\/i>-methylated G22 site and two new sites (Um2 and Am21) were also identified by the irPE method (Supplementary Fig. 6b,c<\/a>). In addition, we identified 34 bktRNA1:U12 chimeras in the FBL CLASH sequencing data (Fig. 3g<\/a> and Supplementary Fig. 7a,b<\/a>). IF and FISH experiments further confirmed that bktRNA1 co-localized with FBL (Supplementary Fig. 8a,b<\/a>). Moreover, knockdown of FBL in HCT116 cells decreased the level of 2\u2032-O<\/i>-methylation of A8 (Am8) in U12 snRNA as well as the level of 2\u2032-O<\/i>-methylation of G22, which is installed by FBL30<\/a><\/sup> (Supplementary Fig. 8c\u2013g<\/a>). These results suggest that bktRNA1 might guide FBL to introduce a conserved 2\u2032-O<\/i>-methylation of the A8 residue in U12 snRNA.<\/p>\n

To investigate bktRNA1 participation in the 2\u2032-O<\/i>-methylation of U12 snRNA, we knocked out bktRNA1 in human HCT116 cells via CRISPR\u2013Cas9. We generated four knockout (KO) clones, including the KO-1 and KO-2 clones with deletion of full-length bktRNA1 and the KO-3 and KO-4 clones with deletion of the guide region and K-turn domains at the 3\u2032 end of bktRNA1 (Fig. 3h<\/a> and Supplementary Fig. 9a\u2013c<\/a>). Quantitative polymerase chain reaction (qPCR) and northern blotting analyses failed to detect accumulation of bktRNA1 in any of the four KO cell lines, and the expression of the U12 snRNA was not affected (Fig. 3i,j<\/a>). Likewise, correct processing and efficient accumulation of CWF19L1 mRNA in bktRNA1-KO cells were verified by RT\u2013PCR and Sanger sequencing (Supplementary Fig. 9d,e<\/a>). The 2\u2032-O<\/i>-methylation state of U12 at A8 was monitored by the irPE method with RNA extracted from wild-type (WT) and KO cells (Fig. 3k<\/a>). Comparing with WT HCT116 cells, the primer extension stop signal was absent from A8 in KO cells, indicating that 2\u2032-O<\/i>-methylation of A8 was completely inhibited in the absence of bktRNA1 (Fig. 3k<\/a>).<\/p>\n

We further restored bktRNA1 expression in KO-4 cells (Supplementary Fig. 9f,g<\/a>) and found that stable expression of full-length bktRNA1, but not partial bktRNAs, in the KO cells (Fig. 3l<\/a>) largely restored the primer extension stop signals at 1\u2009nt downstream of A8 (Fig. 3m<\/a>), indicating that ectopically expressed bktRNA1 fully restored A8 methylation of U12 snRNA. Additionally, bktRNA1 mutants with mutation or truncation of the backward K-turn domain, as well as substitution in the residue that directly interacts with the A8 residue of U12 snRNA, could completely inhibit or significantly diminish A8 methylation of U12 snRNA (Supplementary Figs. 10<\/a> and 11<\/a> and Supplementary Note 2<\/a>). Moreover, artificial guide RNAs that contained the guide sequences of bktRNA1 could induce 2\u2032-O<\/i>-methylation of A8 in U12 snRNA (Supplementary Fig. 12<\/a> and Supplementary Note 2<\/a>). Collectively, these results demonstrate that bktRNA1 is required for efficient 2\u2032-O<\/i>-methylation of U12 snRNA.<\/p>\n

Depletion of bktRNA1 causes dysregulation of U12-type introns<\/h3>\n

Our finding that bktRNA1 is responsible for A8 2\u2032-O<\/i>-methylation in U12 snRNA raised the question of whether bktRNA1 and its modification are required for the splicing of U12-type introns in human cells. To explore this, we performed strand-specific mRNA-seq in all four bktRNA1-KO cell lines (Extended Data Fig. 6a<\/a>). Global analysis of aberrant intron retention revealed substantially greater retention of U12-type introns in bktRNA1\u2212\/\u2212<\/sup> cells than in WT cells (Fig. 4a\u2013d<\/a> and Supplementary Table 5<\/a>). Depletion of bktRNA1 affected more than 75% of U12-type introns (Fig. 4e<\/a>); at least 37% of the retained introns had significant changes (P<\/i>\u2009<\u20090.05; Extended Data Fig. 6b<\/a>). In contrast, depletion of bktRNA1 had a negligible effect on U2-type introns (Extended Data Fig. 6b<\/a>). Moreover, compared to U2-type introns, U12-type introns showed global and considerable increases in inclusion level in all four bktRNA1\u2212\/\u2212<\/sup> cell lines (Fig. 4f\u2013i<\/a>). Furthermore, these retained U12-type introns exhibited striking overlap across the four bktRNA1-KO cell lines (Fig. 4j<\/a>). Intriguingly, differential expression analysis revealed that the aberrant splicing of U12-type introns did not affect the expression levels of mature RNAs, consistent with the findings of a previous study31<\/a><\/sup> (Extended Data Fig. 6c<\/a>). These results suggest that loss of bktRNA1 has a global effect on the splicing of U12-type introns.<\/p>\n

\n
Fig. 4: Dysregulation of U12-type intron splicing in bktRNA1-depleted cells.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n
\n

a<\/b>\u2013d<\/b>, Dot plots displaying the intron retention levels in a representative pairwise analysis of bktRNA1 WT HCT116 cells versus bktRNA1-deficient cells, including KO-1 (a<\/b>), KO-2 (b<\/b>), KO-3 (c<\/b>) and KO-4 (d<\/b>) cells. The red dots represent U12-type introns and the blue dots represent U2-type introns. e<\/b>, Proportion of aberrantly retained U12-type introns in bktRNA1-deficient cells. The red boxes represent retained introns, and the blue boxes represent unretained introns. f<\/b>\u2013i<\/b>, Cumulative fraction of the inclusion level difference between U12-type and U2-type introns in WT and bktRNA1-deficient cells, including KO-1 (f<\/b>), KO-2 (g<\/b>), KO-3 (h<\/b>) and KO-4 (i<\/b>) cells. The P<\/i> values on the cumulative plots of inclusion level differences were calculated using a two-sided Mann\u2013Whitney\u2013Wilcoxon test. j<\/b>, Venn diagram showing the numbers of overlapping retained introns across four bktRNA1-deficient HCT116 cell lines. k<\/b>, Coverage of RNA-seq reads from WT HCT116 and bktRNA1-deficient cells in the corresponding U12-type and U2-type introns of protein-coding genes (RAF1, UBL5 and ATG3). The read counts were normalized using the reads per kilobase per million mapped reads (RPKM) values and are shown using an identical scale for all samples. The intron types are indicated in the bottom panel. l<\/b>, The ratio of spliced to unspliced pre-mRNA for U12-type (RAF1, UBL5 and ATG3) and U2-type (GAPDH) introns in bktRNA1-deficient and WT HCT116 cells was determined by qPCR. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. m<\/b>, RT\u2013PCR analysis of RAF1, UBL5 and ATG3 U12-type intron retention in bktRNA1-deficient and WT HCT116 cells. The schemes of the spliced forms are shown at the right. The boxes indicate exons flanking the minor intron.<\/p>\n

Source data<\/a><\/p>\n<\/div>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n

To validate the aberrant U12-type intron retention identified by RNA-seq, we applied qPCR and RT\u2013PCR to measure splice products in bktRNA1-KO and control cells. As expected, the abundance of all 17 randomly selected U12-type introns was upregulated in bktRNA\u2212\/\u2212<\/sup> cells as compared to WT cells and the abundance of the U2-type intron in GAPDH<\/i> genes (Fig. 4k\u2013m<\/a> and Extended Data Fig. 6d<\/a>). Notably, re-expression of bktRNA1 or artificial guide RNAs (ART7 and ART10) in KO cells restored correct U12-type splicing (Extended Data Fig. 6e,f<\/a>). In contrast, restoration of the H\/ACA domain (SNORA12) of bktRNA1 did not rescue the deficiency in minor splicing (Extended Data Fig. 6g<\/a>). These data strongly suggest that bktRNA1 and Am8 in U12 snRNA are important for the fidelity of U12-type splicing in human cells.<\/p>\n

bktRNA1 modulates the recruitment of ZCRB1 to minor spliceosome<\/h3>\n

To examine whether the loss of bktRNA1 interferes with the assembly of minor spliceosome components, we performed a RIP assay in which RNA\u2013protein complexes were immunoprecipitated using antibodies against U12-associated proteins (ZCRB1, PDCD7, RNPC3 and SNRNP48)32<\/a><\/sup> (Fig. 5a<\/a> and Extended Data Fig. 7a<\/a>). The amount of snRNAs (U11, U12, U4atac, U5 and U6atac) involved in minor spliceosomes in the co-precipitate were then measured by qPCR and northern blotting in WT and bktRNA1 KO cells. Compared to that in WT cells, we found that the ZCRB1-antibody-bound complex in KO cells had obviously lower amounts of U11 and U12 snRNAs (Fig. 5b<\/a> and Extended Data Fig. 7b,c<\/a>). Additionally, these differences were not due to the relative amounts of proteins because no changes in the protein levels of ZCRB1 were observed upon knockout of bktRNA1 (Fig. 5c<\/a>). Overexpression of bktRNA1 in KO cells largely increased the enrichment of the U11 and U12 snRNAs in the ZCRB1 complex (Fig. 5d,e<\/a> and Extended Data Fig. 7d<\/a>), indicating that ectopically expressed bktRNA1 fully restored the interaction of ZCRB1 with U11 and U12. We next performed U12 chromatin immunoprecipitation by RNA purification (ChIRP) assays followed by western blotting to verify the interaction of ZCRB1 and U12 in vivo (Fig. 5f<\/a> and Supplementary Fig. 13a,b<\/a>). In ChIRP analyses, U12 probes retrieved more ZCRB1 in bktRNA1 WT cells than in KO cells (Fig. 5f<\/a>, upper panel). Conversely, re-expression of bktRNA1 in KO cells markedly increased the amounts of ZCRB1 retrieved by U12-specific probes (Fig. 5f<\/a>, lower panel).<\/p>\n

\n
Fig. 5: Depletion of bktRNA1 affects U12 interactions with ZCRB1.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n
\n

a<\/b>, Western blots showing precipitation with an anti-ZCRB1 antibody in WT and bktRNA1-deficient KO-4 cells. b<\/b>, Native RIP was performed in WT and KO-4 cells using an anti-ZCRB1 antibody or normal IgG antibody, after which qPCR was performed with primers recognizing snRNAs involved in minor spliceosomes. The percentage of RIP-enriched snRNAs was calculated relative to the input RNA. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. NS, no significance. c<\/b>, Western blotting analysis of ZCRB1 in HCT116 WT and bktRNA1\u2212\/\u2212<\/sup> clones. GAPDH served as the loading control. d<\/b>, Detection of ZCRB1 in immunoprecipitates from KO-4 and bktRNA1-rescued (bktRNA1) cells by western blotting analysis. e<\/b>, Native RIP was performed in bktRNA1-deficient and bktRNA1-rescued cells using an anti-ZCRB1 antibody or normal IgG antibody followed by qPCR with primers for snRNAs involved in minor spliceosomes. The percentage of RIP-enriched snRNAs was calculated relative to the input RNA. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. f<\/b>, Western blotting analyses of proteins retrieved by the U12-ChIRP experiments. RNPC3 was used as a positive control, and GAPDH was used as a negative control. NC, negative control probes. g<\/b>, Validation of 2\u2032-O<\/i>-methylated U12-specific binding proteins by RNA affinity chromatography using single-stranded RNA probes with 2\u2032-O<\/i>-methylated (red) or unmethylated (green) adenosine. A schematic workflow (left panel) and western blotting analyses (right panel) of RNA affinity chromatography were used to validate U12-interacting proteins (ZCRB1, RNPC3, PDCD7 and SNRNP48). GAPDH was used as a negative control. h<\/b>, RNA EMSA was used to assess the binding of recombinant ZCRB1 and the indicated U12 RNA probes. The dissociation constants (K<\/i>d<\/sub>, nM) of recombinant ZCRB1 proteins with methylated (Methyl-mutU12) or unmethylated (Unmethyl-mutU12) U12 ssRNA probes are indicated in the lower panel. i<\/b>, Proposed model showing the functions and mechanisms of bktRNA1. 65K is RNPC3; 59K is PDCD7; and 48K is SNRNP48.<\/p>\n

Source data<\/a><\/p>\n<\/div>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n

We further performed RNA affinity chromatography followed by immunoblotting to validate the proteins that interact with methylated U12 probes. We discovered that ZCRB1 selectively bound to 2\u2032-O<\/i>-methylated U12 RNA baits compared to the unmethylated control (Fig. 5g<\/a>). Similarly to endogenous proteins, the purified recombinant ZCRB1 proteins (Supplementary Fig. 13c,d<\/a>) also preferentially bound to the 2\u2032-O<\/i>-methylated U12 probe (Methyl-U12, K<\/i>d<\/sub>\u2009=\u2009340.3\u2009\u00b1\u200948.9\u2009nM) over the unmethylated one (Unmethyl-U12, K<\/i>d<\/sub>\u2009=\u20092,293.0\u2009\u00b1\u2009358.4\u2009nM) in vitro (Fig. 5h<\/a>). We further performed an RNA EMSA analysis with mutant U12 probes, which have a change of A to G (UUGA) in the consensus sequence recognized by ZCRB1 (refs. 33<\/a>,34<\/a><\/sup>), and found that recombinant ZCRB1 exhibited higher affinity toward the mutated U12 probe with methylated guanine ribonucleotide (Methyl-mutU12, K<\/i>d<\/sub>\u2009=\u20091,787.0\u2009\u00b1\u2009211.4\u2009nM) than unmethylated probes (Unmethyl-mutU12, K<\/i>d<\/sub>\u2009=\u200913,029.0\u2009\u00b1\u20094,687.0\u2009nM) (Supplementary Fig. 13e<\/a>).<\/p>\n

We then investigated the minor splicing dysregulation and phenotypic changes upon ZCRB1 knockdown in HCT116 cells (Extended Data Fig. 8a<\/a>). The high similarity between the bktRNA1-KO and ZCRB1-knockdown cells for inducing minor intron retention and functional phenotypes further indicates a close connection between bktRNA1 and ZCRB1 (Extended Data Fig. 8b\u2013m<\/a>, Supplementary Table 6<\/a> and Supplementary Note 3<\/a>). Together, these results suggest that the 2\u2032-O<\/i>-methylation at A8 in U12 snRNA guided by bktRNA1 is crucial for the recruitment of ZCRB1 to the U11\u2013U12 di-snRNP complex as well as for the splicing of U12-type introns (Fig. 5i<\/a>).<\/p>\n

bktRNAs regulate the splicing of local introns<\/h3>\n

Given that bktRNAs are universally located within introns and that their binding partner 15.5K can promote the assembly of minor and major spliceosomes7<\/a><\/sup>, we hypothesized that bktRNAs might alternatively be involved in regulating intron processing and RNA splicing. To test this hypothesis, we performed PEN-seq and mRNA-seq in 15.5K-knockdown cells (Supplementary Fig. 14a,b<\/a>). As expected, the expression of bktRNAs was markedly decreased in 15.5K-knockdown cells compared to control cells (Supplementary Fig. 14c,d<\/a> and Supplementary Table 7<\/a>). Notably, intron retention analysis from the mRNA-seq datasets showed that the silencing of 15.5K resulted in the global alteration of intron splicing (Supplementary Fig. 14e,f<\/a> and Supplementary Table 8<\/a>). Moreover, the inclusion level was considerably different between bktRNAs and the remaining introns upon 15.5K knockdown (Supplementary Fig. 14g,h<\/a>).<\/p>\n

We further validated the splicing efficiency of local introns containing bktRNA1 in WT and KO cells by qPCR and found that the depletion of bktRNA1 resulted in downregulation of local intron splicing but did not affect the distant intron splicing (Fig. 6a<\/a>). Moreover, the impaired splicing efficiency of local introns could not be rescued by re-expression of WT bktRNA1 in either bktRNA1 KO-1 or KO-4 cells (Fig. 6b<\/a>), implying that bktRNA1 might regulate the local splicing of introns by sequence elements of bktRNA1 but not the specific bktRNA1 transcript, which may be similar to some long non-coding RNAs (lncRNAs)35<\/a><\/sup>.<\/p>\n

\n
Fig. 6: Local regulation of intron splicing by bktRNAs.<\/b><\/figcaption>
\n
\"Science<\/picture><\/a><\/div>\n

a<\/b>, The ratio of spliced to unspliced introns containing bktRNA1 (local intron) and another intron from the same host gene (distant intron in CWF19L1) was determined by qPCR in bktRNA1-deficient and WT HCT116 cells. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. b<\/b>, Ratio of spliced to unspliced intron (local) in bktRNA1-rescued HCT116 KO-1 (right panel) and KO-4 (left panel) cells. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. c<\/b>, Schematic of the GFP splicing reporter design. bktRNA cloning and primers were used for splicing efficiency detection. We used the primer pairs sF1\/R1 and usF1\/R1 to detect mature and precursor GFP, respectively. The analysis of puromycin with primer pairs qF1\/qR1 acted as an internal control. d<\/b>, qPCR analysis of the ratio of spliced to unspliced GFP RNA in HEK293T cells transfected with WT and mutant bktRNAs (mutated from CTGA to CTAG). Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. e<\/b>, The ratio of spliced to unspliced intron containing bktRNA161 (local intron) and another intron from the same host gene (distant intron in TEAD1) was determined by qPCR in bktRNA161-mutant and WT HEK293T cells. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20094 biological replicates), two-tailed, paired t<\/i>-test. f<\/b>, Quantification of introns containing bktRNA161 by qPCR in 15.5K immunoprecipitates from WT and mutant bktRNA161 HEK293T cells. U4atac and MT-CO1 were kept as positive and negative controls, respectively. Data are presented as mean values \u00b1 s.e.m. (n<\/i>\u2009=\u20093 biological replicates), two-tailed, paired t<\/i>-test. g<\/b>, Model of local regulation of intron splicing by bktRNAs. The interaction between the backward K-turn structure of bktRNAs and 15.5K negative (for example, bktRNA161) or positive (for example, bktRNA1) regulates the splicing of the local intron containing bktRNAs. In contrast, disruption of the backward K-turn structure in bktRNA161 or bktRNA1 relieves the local regulation of intron splicing.<\/p>\n<\/div>\n

Full size image<\/span><\/a><\/p>\n<\/figure>\n<\/div>\n

To test the role of bktRNAs in the local regulation of intron splicing, we constructed GFP reporters that harbor bktRNAs with WT or backward K-turn motif mutants (mutated from CUGA to CUAG, Mut) (Fig. 6c<\/a>). Consistent with the local regulation of intron splicing by bktRNA1, WT bktRNA1 in the GFP reporter intron also promoted the splicing of local introns, whereas mutation of the 5\u2032 motif (bktRNA1-MM) or deletion of the backward K-turn motif significantly decreased the splicing of local introns (Fig. 6d<\/a>). By applying this GFP reporter to ten randomly selected bktRNAs, we found that eight bktRNAs can regulate the splicing of local introns (Fig. 6d<\/a>). For example, WT bktRNA168 enhanced the local splicing in comparison to the mutant type (Fig. 6d<\/a>), whereas bktRNA2 with a backward K-turn-mutant increased the splicing of local introns (Fig. 6d<\/a>). These results suggest that the backward K-turn motif of bktRNAs is important for the splicing of local introns.<\/p>\n

To confirm the local regulation of bktRNAs in intron splicing in vivo, we mutated the 5\u2032 motif (from CTGA to CTAG) of the backward K-turn structure in bktRNA161<\/i> with prime editing36<\/a><\/sup> (Supplementary Fig. 15a<\/a>). We generated seven homozygous mutants (from mut161-1 to mut161-7) based on Sanger sequencing (Supplementary Fig. 15b<\/a>). In line with the GFP reporter assays (Fig. 6d<\/a>), the splicing of the local intron containing bktRNA161, compared to the distant intron, was enhanced upon mutation of the backward K-turn structural motif of bktRNA161 in all mutants (Fig. 6e<\/a>). Moreover, the binding affinity of the bktRNA161-containing intron with the 15.5K protein was substantially reduced in bktRNA161 mutant HEK293T cells (Fig. 6f<\/a> and Supplementary Fig. 15c<\/a>). Together, these results strongly support that bktRNAs are involved in the local regulation of intron splicing in a backward K-turn motif-dependent manner (Fig. 6g<\/a>).<\/p>\n<\/div>\n<\/div>\n

\n

Discussion<\/h2>\n
\n

In this study, we discovered a class of bktRNAs that act as binding partners for the K-turn binding protein 15.5K by developing the RIP-PEN-seq method, and we revealed that bktRNA1 can act as a post-transcriptional regulator of minor splicing and that bktRNAs can regulate the splicing of local introns. Although U12-type introns constitute ~0.5% of the introns in the human genome, they exist in approximately 700\u2013800 genes involved in various cellular processes37<\/a>,38<\/a><\/sup>. Mutations in genes encoding the components of the U12 spliceosome, including both snRNAs and interacting proteins, have been shown to cause a series of developmental disorders in humans39<\/a>,40<\/a>,41<\/a>,42<\/a>,43<\/a><\/sup>. For example, mutation in U12 is associated with cerebellar ataxia44<\/a><\/sup>. These findings raise an interesting question: whether aberrant expression of bktRNA1 could lead to some of the abovementioned developmental disorders in humans. Coincidentally, mutation of CWF19L1<\/i>, the host gene of bktRNA1<\/i>, has been considered as cause of autosomal recessive cerebellar ataxia (ARCA)45<\/a>,46<\/a>,47<\/a><\/sup>, consistent with the phenotype caused by the mutation in U12 snRNAs44<\/a><\/sup>. Given that the loss of bktRNA1 leads to global splicing dysregulation of U12-type introns, we hypothesized that mutation of CWF19L1<\/i> might affect the expression and processing of bktRNA1 to cause ARCA disease. Therefore, future works might address whether bktRNA1<\/i> is associated with ARCA.<\/p>\n

Although backward K-turn structural motifs have been identified in few ncRNAs (for example, 23S rRNA from Archaea<\/i>)2<\/a><\/sup>, these motifs are not composed of consensus sequences and are located within different positions of ncRNAs. In contrast to these few ncRNAs2<\/a><\/sup>, our bktRNAs identified from mammals all have two consensus sequence elements located within backward K-turn motifs. Importantly, these two consensus sequence elements are located at fixed positions in bktRNAs (Fig. 1c\u2013e<\/a>). Moreover, we found that these backward K-turn structural motifs of bktRNAs are indispensable for avoiding degradation by exonucleases, the processing and maturation of bktRNAs (Supplementary Fig. 10<\/a> and Supplementary Note 2<\/a>) and the local regulation of intron splicing (Fig. 6<\/a>).<\/p>\n

We used the bktRNA1 gene as a model to explore the regulatory function of bktRNAs (Fig. 5i<\/a>). Moreover, the consensus structural motifs of bktRNAs and their dynamic, cell-specific\/tissue-specific expression patterns have shed light on their roles in gene expression regulation. We found that bktRNA1 and other bktRNAs are involved in the local regulation of intron splicing (Fig. 6g<\/a>). Notably, the backward K-turn structural motifs are indispensable for the local regulation of intron splicing by bktRNAs (Fig. 6<\/a>). Interestingly, bktRNAs can act not only as positive regulators but also as negative regulators of the local splicing of introns (Fig. 6g<\/a>), suggesting that bktRNAs may interact with different factors to regulate the splicing of local introns. Thus, it would be very interesting to explore any cofactors (for example, splicing factors) that control the local regulation of intron splicing by these bktRNAs.<\/p>\n

In summary, our data strongly suggest that bktRNAs form a class of post-transcriptional regulators that probably help RBPs recognize substrates through base-pairing interactions (for example, bktRNA1) and generally function in modulating their binding partner 15.5K to regulate the splicing of local introns. Our studies reveal that a single methylation alteration caused by the loss of bktRNA1 contributes to global splicing dysregulation and that bktRNAs can act as local regulators in the splicing of introns. Overall, these discoveries add another layer of complexity to the control of gene expression that involves crosstalk among RNA splicing, RNA methylation and bktRNAs.<\/p>\n<\/div>\n<\/div>\n

\n

Methods<\/h2>\n
\n

Construction of RIP-PEN-seq and PEN-seq libraries<\/h3>\n

For conventional RIP-seq, the immunoprecipitated RNAs (especially for RNAs with length \u226550\u2009nt) were often fragmented and then subjected to RNA-seq library construction with random primer-based reverse transcription49<\/a><\/sup>. As a result, the conventional RIP-seq method cannot identify the full length of immunoprecipitated RNA, and, thus, conventional RIP-seq cannot discover the precise positions of motifs at RNAs. In this study, we used dual RNA adapters to capture both ends of ncRNAs. This strategy and the following strategies allowed us to identify the full-length sequences of ncRNAs (from 20\u2009nt to 500\u2009nt) and facilitated the discovery of bktRNAs with characteristic motifs and structures presented at RNA ends.<\/p>\n

RIP-PEN-seq and PEN-seq (Fig. 1b<\/a> and Extended Data Fig. 1a<\/a>) employed several strategies to achieve high-performance ncRNA full-length sequencing. First, to enable the high ligation efficiency and quantitative full-length library, specific randomized barcode-containing (4N, N represents for A, T or U, G, C) 5\u2032 adapters and 3\u2032 adapters were synthesized and directly ligated to RNAs. Second, to increase the diversity of RNA species and discover low-expression RNAs, we used RNase H or Ribo-Zero strategies to remove high-abundance rRNAs and\/or snRNAs50<\/a>,51<\/a><\/sup>. It is worth noting that the Ribo-Zero kit from Illumina had been discontinued, and the nuclease-mediated depletion displays huge biases in sRNA-seq following the traditional procedure52<\/a><\/sup>. Here, we performed rRNA\/snRNA depletion after adaptors ligation and used single-stranded DNA (ssDNA) exonuclease RecJf rather than DNase I to remove the probes. Third, to avoid the widespread mispriming artifacts53<\/a><\/sup> and obtain full-length cDNAs, we used a nested reverse transcription primer during reverse transcription, along with a full-length reverse PCR primer for amplification. Finally, to overcome the obstacles of RNA modifications and stable RNA secondary structures in cDNA synthesis, we used SuperScript IV Reverse Transcriptase, which has a high level of processivity and thermostability, to generate cDNAs.<\/p>\n

RIP experiments were performed using the protocol outlined in a previous report54<\/a>,55<\/a><\/sup>. In brief, for each immunoprecipitate, 1\u2009\u00d7\u2009107<\/sup> cells stably expressing FLAG-15.5K were pelleted; resuspended in an equal volume of ice-cold polysomal lysis buffer (10\u2009mM HEPES pH 7.0, 100\u2009mM KCl, 5\u2009mM MgCl2<\/sub>, 0.5% NP-40) supplemented with 1\u2009mM DTT, 100\u2009U\u2009ml\u22121<\/sup> RNase inhibitor (Takara), 1\u00d7 protease inhibitor cocktail (Roche) and 0.4\u2009mM RVC (New England Biolabs (NEB)); incubated on ice for 15\u2009minutes; snap frozen in liquid nitrogen; and stored at \u221280\u2009\u00b0C. After the lysate was thawed on ice and centrifuged at 15,000g<\/i> for 15\u2009minutes, the supernatant was pre-cleared with Dynabeads Protein G (Invitrogen) at 4\u2009\u00b0C for 1\u2009hour. Ten microliters of the supernatant was saved as input and mixed with 1\u2009ml of RNAzol. Subsequently, a 100-\u03bcl aliquot of the cell extract was diluted with 900\u2009\u03bcl of NT2 buffer (50\u2009mM Tris pH 7.4, 150\u2009mM NaCl, 1\u2009mM MgCl2<\/sub>, 0.05% NP-40) supplemented with 1\u2009mM DTT, 100\u2009U\u2009ml\u22121<\/sup> RNase inhibitor (Takara), 1\u00d7 protease inhibitor cocktail (Roche) and 20\u2009mM EDTA. The extracts were mixed with 5\u2009\u03bcg of an anti-FLAG antibody (Proteintech) or a normal IgG antibody (Millipore) and incubated overnight at 4\u2009\u00b0C with end-over-end rotation. The next day, the RNP\/antibody complex was precipitated by incubation with Dynabeads Protein G at 4\u2009\u00b0C for 3\u2009hours. The beads were collected and washed five times in NT2 buffer. After the final wash, the beads were resuspended in 1\u2009ml of RNAzol and saved as RIP RNA or IgG RNA.<\/p>\n

For PEN-seq, total RNA was isolated from cell lines with RNAzol, treated with RQ1 DNase (Promega) and purified by RNA Clean & Concentrator-5 (Zymo Research). The samples were then ligated to 3\u2032 randomized RNA adaptor with T4 RNA ligase 2 truncated KQ (NEB) in 1\u00d7 T4 RNA ligase reaction buffer supplemented with 12.5% PEG 8000 at 16\u2009\u00b0C for 18\u2009hours. The excessive adaptors were digested with 100\u2009U of 5\u2032 deadenylase (NEB) at 30\u2009\u00b0C for 1\u2009hour, incubated with 2\u2009\u03bcg of Escherichia coli<\/i> single-strand DNA-binding protein (Promega) on ice for 30\u2009minutes and subjected to ssDNA digestion with 60\u2009U of RecJf (NEB) at 37\u2009\u00b0C for another hour. The ligated RNA was ligated to a 5\u2032 randomized RNA adaptor using T4 RNA Ligase 1 (NEB) in 1\u00d7 T4 RNA ligase reaction buffer supplemented with 1\u2009mM ATP at 16\u2009\u00b0C for 18\u2009hours. The ligated RNA was column-purified by RNA Clean & Concentrator-5 (Zymo Research) and subjected to RNase H-based ribosomal RNA depletion according to a standard protocol with some modifications51<\/a><\/sup>. For rRNA\/snRNA probe digestion, DNase I was replaced by RecJf exonuclease to prevent digestion of the 3\u2032 randomized RNA adaptor. After column-based purification, the rRNA\/snRNA-depleted RNA was reverse transcribed using SuperScript IV Reverse Transcriptase (Invitrogen) in 1\u00d7 reverse transcription buffer with truncated reverse transcription primer to eliminate mispriming artifacts. cDNA was amplified using Phusion High-Fidelity DNA Polymerase (Thermo Fisher Scientific) in 1\u00d7 Phusion HF buffer (Thermo Fisher Scientific) with primers (the RP1 forward primer and RPI-X reverse primer). The PCR-amplified cDNA was gel purified using a 4% NuSieve 3-1 Agarose gel (Lonza) to remove adapter dimers and sequenced on an Illumina HiSeq X Ten platform with paired-end reads (150\u2009bp or 125\u2009bp) at Annoroad Gene Technology Company. The adapters and primers are listed in Supplementary Table 9<\/a>.<\/p>\n

RIP-PEN-SHAPE-MaP library construction<\/h3>\n

RIP-PEN-SHAPE-MaP was performed on HEK293T-F-15.5K cells as previously described21<\/a><\/sup> with some modifications. In brief, RIP experiments were carried out as described above until the last NT2 buffer wash. After the last wash, the beads were incubated with SHAPE modification buffer (100\u2009mM HEPES pH 8.0, 6\u2009mM MgCl2<\/sub>, 150\u2009mM NaCl, 50\u2009mM NAI-N3<\/sub> (MedChemExpress, HY-103006))56<\/a><\/sup> at 37 \u00b0C for 12\u2009minutes on a thermomixer at 1,000\u2009r.p.m. The control group was prepared using DMSO following the same procedure as described above. RNA was extracted using RNAzol, treated with DNase I and purified with RNA Clean & Concentrator-5 (Zymo Research) in accordance with the manufacturer\u2019s protocol. Then, approximately 300\u2009ng of NAI-N3-treated or DMSO-treated RNA sample was subjected to the PEN-seq procedure with some modifications. (1) The 3\u2032 randomized RNA adaptor and 5\u2032 randomized RNA adaptor were replaced with randomized adaptors containing six degenerate nucleotides. (2) To increase the read coverage on bktRNAs, we designed probes57<\/a><\/sup> for RNase H-based high-abundance snoRNA subtraction in addition to rRNA and snRNA depletion. The high-abundance snoRNAs were determined according to the top 1,000 reads in 15.5K RIP-PEN-seq datasets. (3) Reverse transcription was performed using SuperScript II Reverse Transcriptase (Thermo Fisher Scientific) in MnCl2<\/sub>-containing 1\u00d7 reverse transcription buffer (50\u2009mM Tris-HCl pH 8.0, 75\u2009mM KCl, 6\u2009mM MnCl2<\/sub>, 10\u2009mM DTT)20<\/a><\/sup>. The whole reaction system in the tube was incubated at 42\u2009\u00b0C for 3\u2009hours. (4) For cDNA purification, Exo I (NEB) was used to digest excess reverse transcription primers at 37 \u00b0C for 15\u2009minutes. Then, 7\u2009\u03bcl of 1\u2009M NaOH and 5\u2009\u03bcl of 0.5\u2009M EDTA per 20-\u03bcl reaction volume were added to remove RNA templates at 70\u2009\u00b0C for 12\u2009minutes, after which cleanup of the cDNA was performed with Oligo Clean & Concentrator (Zymo Research). (5) Library amplification was performed using NEBNext Ultra II Q5 Master Mix (NEB) and was programmed as follows: stage I: 98\u2009\u00b0C for 30\u2009s; stage II: 98\u2009\u00b0C for 10\u2009s, 65\u2009\u00b0C for 75\u2009s, with a limited cycle; stage III: 65\u2009\u00b0C for 5\u2009minutes. The number of cycles was 16\u201320. The adapters, primers and probes are listed in Supplementary Table 9<\/a>.<\/p>\n

Mapping of 2\u2032-O<\/i>-methylation residues<\/h3>\n

The positions of 2\u2032-O<\/i>-methyls were identified by dNTP concentration-dependent primer extension (PE) assays. The primer used for U12 primer extension (irPE-U12) was synthesized by Sangon. To label the azide-modified oligonucleotides with IR dye, 0.4\u2009nmol of oligonucleotides was mixed with 10\u2009nmol IRDye 680RD DBCO (LI-COR Biosciences) in 1\u00d7 PBS (Thermo Fisher Scientific) at 25\u2009\u00b0C overnight. The IR dye-labeled oligonucleotides were then purified with Oligo Clean & Concentrator kits (Zymo Research). The 5\u2032 end labeled primers were extended by SuperScript III Reverse Transcriptase (Invitrogen) in the presence of decreasing amounts (0.5\u2009mM, 0.05\u2009mM and 0.01\u2009mM) of dNTPs using 30\u2009\u03bcg of total RNA to determine modifications in U12 snRNA. Sequence analysis was carried out according to the dideoxy method with a Sequenase Version 2.0 DNA Sequencing Kit (USB). The products were separated on a 10% urea-polyacrylamide gel and visualized on a LI-COR Odyssey Infrared Imager with Odyssey software version 3.0.<\/p>\n

ChIRP<\/h3>\n

ChIRP was performed according to previous studies58<\/a>,59<\/a>,60<\/a><\/sup> with some modifications. HCT116 cells were crosslinked with 3% formaldehyde at room temperature for 30\u2009minutes, followed by 125\u2009mM glycine quenching at room temperature for 5\u2009minutes. Cells were centrifuged, and the pellet was washed with ice-cold PBS and then snap frozen in liquid nitrogen and stored at \u221280\u2009\u00b0C. Crosslinked cells (~2\u2009\u00d7\u2009107<\/sup>) were resuspended with 1\u2009ml of ChIRP pre-digestion buffer (20\u2009mM Tris-HCl pH 7.5, 5\u2009mM MgCl2<\/sub>, 0.5\u2009mM CaCl2<\/sub>, 0.5% Triton X-100) supplemented with 1\u00d7 protease inhibitor cocktail (Roche), 1\u2009mM DTT, 100\u2009U\u2009ml\u22121<\/sup> RNase inhibitor, 10\u2009mM RVC (NEB)n and 12\u2009U of DNase I (Promega) and rotated at 37\u2009\u00b0C for 10\u2009minutes. The reaction was stopped with 20\u2009mM EDTA, and cells were pelleted, washed once with 1\u00d7 PBS and resuspended in ChIRP cell lysis buffer (50\u2009mM Tris-HCl pH 7.0, 10\u2009mM EDTA, 1% SDS) with 1\u00d7 protease inhibitor cocktail (Roche), 1\u2009mM DTT and 100\u2009U\u2009ml\u22121<\/sup> RNase inhibitor. Cells were sonicated using Sonic Vibra-Cell VCX130 for 20 cycles with the following conditions: 30% amplitude, 10\u2009s on, 10\u2009s off. After centrifugation at 16,000g<\/i> for 10\u2009minutes at 4\u2009\u00b0C, the supernatant was collected and pre-cleared with Dynabeads MyOne Streptavidin C1 beads (Thermo Fisher Scientific) for 30\u2009minutes at 37\u2009\u00b0C. For hybridization, the pre-cleared cell lysates were diluted in two times volume of ChIRP hybridization buffer (50\u2009mM Tris-HCl, pH 7.0, 750\u2009mM NaCl, 1\u2009mM EDTA, 1% SDS, 15% formamide) with 1\u00d7 protease inhibitor cocktail (Roche), 1\u2009mM DTT and 100\u2009U\u2009ml\u22121<\/sup> RNase inhibitor and incubated with 100\u2009pmol of U12-targeted biotinylated probe mix (ChIRP-U12-1, 2, 3) or negative control probe (ChIRP-NC) at 37\u2009\u00b0C overnight. The next day, 100\u2009\u03bcl of pre-washed Dynabeads MyOne Streptavidin C1 beads were added and rotated for an additional 30\u2009minutes. The beads were washed five times with ChIRP washing buffer (2\u00d7 SSC, 0.5% SDS) at 37\u2009\u00b0C. A quarter of the ChIRP products were used for RNA isolation, and the remaining were subjected to protein elution. For RNA isolation, input and ChIRP samples were resuspended in ChIRP RNA-isolation buffer (10\u2009mM Tris-HCl pH 7.0, 100\u2009mM NaCl, 1\u2009mM EDTA, 0.5% SDS) supplemented with 1,000\u2009mg\u2009ml\u22121<\/sup> proteinase K and incubated at 50\u2009\u00b0C for 45\u2009minutes, followed by 95\u2009\u00b0C for 15\u2009minutes and RNAzol extraction. Eluted RNA was subjected to ChIRP-qPCR for detection of related transcripts. For protein extraction, input and ChIRP samples were resuspended in ChIRP biotin elution buffer (7.5\u2009mM HEPES pH 7.5, 12.5\u2009mM D-biotin, 1.5\u2009mM EDTA, 75\u2009mM NaCl, 0.15% SDS, 0.075% Sarkosyl, 0.02% Na-deoxycholate), mixed at 25\u2009\u00b0C for 20\u2009minutes and at 65\u2009\u00b0C for 10\u2009minutes. After collecting the eluents, the beads were eluted again. The two eluents were pooled and precipitated with 20\u2009\u03bcg\u2009ml\u22121<\/sup> protein carrier (BSA) and four times volume of ice-cold acetone at \u221220\u2009\u00b0C overnight. The next day, proteins were pelleted and dissolved in 1\u00d7 SDS loading buffer (Beyotime) and boiled at 95\u2009\u00b0C for 30\u2009minutes. Final protein samples were subjected to SDS-PAGE for western blotting. The probes used for ChIRP are listed in Supplementary Table 9<\/a>.<\/p>\n

Expression and purification of ZCRB1 protein<\/h3>\n

We use the E. coli<\/i> system to produce and purify specific ZCRB1 protein following the HaloTag protein purification system (Promega) with pieces of modifications. In brief, pH6HTN-ZCRB1 vectors that expressed HaloTag fusion ZCRB1 protein were transformed into E. coli Transetta DE3<\/i> chemically competent cells (TransGen Biotech). A fresh colony was inoculated into 2YT medium containing the antibiotic and grown at 37\u2009\u00b0C for overnight. On the second day, the starter culture was diluted as 1:100 in fresh 2YT medium with antibiotic and cultured at 37\u2009\u00b0C for about 6\u2009hours, reaching the logarithmic growth phase. Then, 1\u2009ml of culture was seeded into 100\u2009ml of 2YT medium with IPTG at 0.4\u2009nM final concentration and grown for 4\u2009hours at 25\u2009\u00b0C for harvesting. According to the HaloTag protein purification system instructions, cells in 100\u2009ml of medium were harvested by centrifugation at 4,000g<\/i> for 10\u2009minutes at 4\u2009\u00b0C, and then the cell pellets were resuspended and lysed by sonication (Sonic Vibra-Cell VCX130) in 10\u2009ml of Halo purification buffer (50\u2009mM HEPES pH 7.5, 150\u2009mM NaCl, 0.5\u2009mM EDTA, 0.005% NP-40) plus with 1\u2009mM DTT, complete protease inhibitor cocktail (Roche), 1\u2009mg\u2009ml\u22121<\/sup> lysozyme and 0.01\u2009U of DNase I for 5\u2009s on, 5\u2009s off, for 36 cycles with 35% amplitude. After centrifugation at 10,000g<\/i> for 15\u2009minutes at 4\u2009\u00b0C, the supernatant was transferred for HaloTag protein purification by HaloTag Resin, which had been pre-washed using HaloTag purification buffer. With binding the HaloTag-ZCRB1-fused protein and centrifugation, the supernatant (the sample flowthrough) was discarded. By washing with HaloTag purification buffer for a total of three times, the HaloTag Resin bound with ZCRB1 was clove with TEV protease cleavage solution (60\u2009\u03bcl of TEV protease supplied in 1\u2009ml of HaloTag purification buffer), which was removed from ZCRB1 elution by HisLink Resin. Finally, 2\u2009ml of eluted ZCRB1 was concentrated with an Amicon Ultra-4 Centrifugal Filter device (10\u2009K<\/i>d<\/i><\/sub>) for 30\u2009minutes up to 100-\u03bcl volume. The final ZCRB1 protein and the reservations at every purified step were analyzed by Coomassie blue staining and western blot, simultaneously.<\/p>\n

RNA electrophoretic mobility shift assays<\/h3>\n

RNA oligonucleotides labeled Cy5 were synthesized by RiboBio. The REMSA was carried out per the instructions of the LightShift Chemiluminescent RNA EMSA Kit (Thermo Fisher Scientific) with some modifications. In brief, RNA probes in annealing buffer (20\u2009mM Tris pH 7.5, 200\u2009mM KCl) were denatured and annealed by incubation at 95\u2009\u00b0C for 4\u2009minutes and then slowly cooled down to 25\u2009\u00b0C at 0.05\u2009\u00b0C per second. Recombinant ZCRB1 was diluted to concentration series of 0\u2009nM, 50\u2009nM, 100\u2009nM, 200\u2009nM and 300\u2009nM in Halo purification buffer. For each reaction, 1 \u03bcl of RNA probes with 6.25\u2009nM final concentration and 1\u2009\u03bcl of ZCRB1 protein (0\u2009nM, 50\u2009nM, 100\u2009nM, 200\u2009nM and 300\u2009nM) were incubated in 2\u2009\u03bcl of 10\u00d7 REMSA Binding Buffer (supplied in the RNA EMSA Kit) plus 5% glycerol and 2\u2009\u03bcg of transfer RNA (tRNA) at room temperature for 30\u2009minutes after adding 5\u2009\u03bcl of 5\u00d7 loading buffer (50\u2009mM HEPES pH 7.5, 80% glycerol, 0.25% bromophenol blue) after separation by native PAGE. The fluorescence signal was visualized by Odyssey Imaging Systems and quantified by Image Studio, and then the dissociation constant (K<\/i>d<\/sub>) was calculated with nonlinear curve fitting (function one-site specific binding) using GraphPad Prism with Y\u2009=\u2009Bmax<\/sub>\u2009\u00d7\u2009X\u2009\/\u2009(K<\/i>d<\/sub>\u2009+\u2009X), where Y is the ratio of [RNA\u2013protein]\u2009\/\u2009([free RNA]\u2009+\u2009[RNA\u2013protein]), X is the input protein concentration and Bmax<\/sub> is set to 1.<\/p>\n

Identification of ktRNAs from RIP-PEN-seq or PEN-seq data<\/h3>\n

We clustered the overlapping paired-end reads that mapped to the genomes. Specifically, we defined the start and end sites with the highest abundance in the cluster as the transcription start site (TSS) and transcription terminal site (TTS), respectively. The sequences between TSSs and TTSs co-occurring within 500\u2009bp were extracted as candidate transcripts. Each candidate transcript with a precise TSS and TTS was input to kturnSeeker. A new software program, kturnSeeker, was developed to discover RNAs with fktRNAs and bktRNAs. kturnSeeker searches A\u2022G, G\u2022A and G\u2022C base pairs (the NC-stem), the 3-nt bulge and a section with regular base pairing (the C-stem) step by step in candidate sequences and scores the C-stem with a strict scoring scheme (a G\u2022C pair is bonus 3, an A\u2022U pair is bonus 2 and mismatch is minus 5). Only those with a total score (the C-stem) greater than 7 were retained as candidate ktRNAs. kturnSeeker can, for the first time, screen bktRNAs as well as fktRNAs by reverse searching the K-turn structure. To identify novel ktRNAs, ktRNAs overlapping with known annotations were discarded. Moreover, all novel ktRNAs must be present in at least two sequencing libraries.<\/p>\n

Annotation of ktRNAs<\/h3>\n

The genome sequence of humans (hg38) was downloaded from the UCSC Genome Browser site61<\/a><\/sup>. Human gene annotations were acquired from GENCODE62<\/a><\/sup> in August 2018. The repeat elements in RepeatMasker were downloaded from the UCSC Genome Browser site61<\/a><\/sup>. The sequences and annotation data for canonical human and mouse box C\/D snoRNAs were downloaded from snoRNA-LBME-db63<\/a><\/sup>, deepBase64<\/a><\/sup>, GENCODE62<\/a><\/sup>, snoRNAome65<\/a>,66<\/a><\/sup> and RefSeq67<\/a><\/sup> in August 2018. All C\/D RNAs were intersected with canonical C\/D box snoRNAs using BEDTools software68<\/a><\/sup>, and only C\/D box RNAs that did not overlap with canonical C\/D box snoRNAs were identified as novel ktRNAs. All novel ktRNAs were further intersected with canonical gene annotations using BEDTools software68<\/a><\/sup>.<\/p>\n

Identification of homologous bktRNA families in other species with different evolutionary clades<\/h3>\n

We constructed homologous bktRNA families based on DNA sequence similarity. To identify evolutionary conservation of human bktRNA families, we used the UCSC liftOver tool (version 396)61<\/a><\/sup> and BLAST+ (version 2.11.0)69<\/a><\/sup> to obtain the homologous regions and sequences of the human bktRNAs in other species with different evolutionary clades. The homologous sequences from other species were input into kturnSeeker software to identify candidate bktRNAs. The genomes included the Pan troglodytes<\/i> reference genome (chimpanzee, panTro6), the Pan paniscus<\/i> reference genome (bonobo, panPan3), the Gorilla gorilla gorilla<\/i> reference genome (gorilla, gorGor6), the Pongo pygmaeus abelii<\/i> reference genome (orangutan, ponAbe3), the Macaca mulatta<\/i> reference genome (rhesus, rheMac10), the Mus musculus<\/i> reference genome (mouse, mm10), the Monodelphis domestica<\/i> reference genome (opossum, monDom5), the Gallus gallus<\/i> reference genome (chicken, galGal6) and the Danio rerio<\/i> reference genome (zebrafish, danRer10).<\/p>\n

Processing of 15.5K and FBL CLASH sequencing data<\/h3>\n

First, Cutadapt (version 2.8)70<\/a><\/sup> was used to cut the adapters of the paired-end reads. The paired-end reads without adapters were further merged using FLASH software (version 1.2.11)71<\/a><\/sup> with the following parameters: -z -m 10 -M 135 -x 0.25. fastx_collapser software was used to collapse the reads with the default parameters for removal of PCR duplicates. Reads with fully identical sequences, including the same barcode and insert fragment, were marked as PCR duplicates, and only one of the identical sequences was retained for subsequent analysis. Reads shorter than 30\u2009nt were removed.<\/p>\n

Differential splicing analysis<\/h3>\n

Adaptor sequences were trimmed from raw RNA-seq data using Cutadapt (version 2.8)70<\/a><\/sup>. The clean reads were mapped to the reference genome (hg38) using STAR software (2.7.1a)72<\/a><\/sup> with the genome index built from GENCODE version 32 (ref. 73<\/a><\/sup>) annotation and with the following additional parameters: \u2013alignEndsType EndToEnd\u2013outSAMstrandField intronMotif\u2013outFilterMismatchNmax 5. The genomic coordinates of introns were extracted from GENCODE version 32 (ref. 73<\/a><\/sup>) annotation, and duplicates were further removed. For each intron, the numbers of reads that supported either the spliced isoform or the retained isoform were counted. Reads that skipped the intron and spanned at least 10\u2009bp in both exons were thought to support the spliced isoform, and reads that included the intron with at least 10\u2009bp overlapping between the exon side and the intron side were thought to support the retained isoform. The read counts were normalized by effective length, defined as the number of possible positions for supporting reads, to estimate the abundance of the spliced isoform and the retained isoform. Then, the inclusion level was calculated as the abundance of the retained isoform divided by the abundance of both isoforms. The statistical method and codes from rMATS74<\/a><\/sup> were used to look for introns with significant inclusion level differences between WT and KO cells. For each pair of compared groups, introns with a sum of the number of reads supporting spliced isoforms or retained isoforms less than 20 in either sample were filtered before statistical testing to remove potential false positives. The annotation of minor introns was downloaded from MIDB38<\/a><\/sup>. The P<\/i> value on the cumulative plots of inclusion level differences were calculated using two-sided Mann\u2013Whitney\u2013Wilcoxon tests.<\/p>\n

Statistics and reproducibility<\/h3>\n

Data are presented as the mean values \u00b1 s.e.m. We used paired Student\u2019s t<\/i>-tests for comparisons between the two experimental groups. All statistics were performed using GraphPad Prism 6 or R (3.6.3). The numbers of biological replicates for the experiment are indicated in the figure legends. At least four independent experiments of RIP-PEN-seq were performed for both HEK293T-FLAG-15.5K and Hepa1-6-FLAG-15.5K cells. Four independent experiments of RIP-PEN-SHAPE-MaP were performed for HEK293T-FLAG-15.5K cells. The PEN-seq analysis for cellular fractions in HEK293T and HCT116 cells was performed four times, respectively. Two or three independent sets of HEK293T, HCT116, U-87 MG, HepG2, HeLa and K562 RNA samples were used for PEN-seq analysis. All images of polyacrylamide gels, agarose gels and northern blots are representative of at least two biological replicates. No statistical method was used to predetermine the sample size.<\/p>\n

Reporting summary<\/h3>\n

Further information on research design is available in the Nature Portfolio Reporting Summary<\/a> linked to this article.<\/p>\n<\/div>\n<\/div><\/div>\n

\n
\n

Data availability<\/h2>\n

All sequencing data that support the findings of this study have been deposited in the National Center for Biotechnology Information\u2019s Gene Expression Omnibus with the following accession numbers: GSE160970<\/a> for all HEK293T 15.5K RIP-PEN-seq; GSE182757<\/a> for all Hepa1-6 RIP-PEN-seq; GSE160636<\/a> for FBL and 15.5K CLASH-seq in HEK293T cells; GSE160887<\/a> for PEN-seq in HCT116, U-87 MG, Hela, HEK293T, HepG2 and K562 cells; GSE186849<\/a> for PEN-seq in 15.5K knockdown HEK293T cells; GSE182843<\/a> for PEN-seq in HEK293T and HCT116 cell fractions; GSE160515<\/a> for RNA-seq in bktRNA1 KO HCT116 cells; GSE182830<\/a> for RNA-seq in ZCRB1 knockdown HCT116 cells; GSE182759<\/a> for RNA-seq in 15.5K knockdown HEK293T cells; and GSE220470<\/a> for RIP-PEN-SHAPE-MaP in HEK293T cells. All data are available in the manuscript and in Supplementary Information<\/a> and Source data<\/a> files. There are no restriction on data availability. Source data<\/a> are provided with this paper.<\/p>\n<\/div>\n

\n

Code availability<\/h2>\n

The program kturnSeeker was written in the C++ programming language and is available from GitHub with no restrictions or conditions on access: https:\/\/github.com\/sysu-software\/kturnSeeker<\/a>.<\/p>\n<\/div>\n

\n

References<\/h2>\n
\n
    \n
  1. \n

    Butcher, S. E. & Pyle, A. M. The molecular interactions that stabilize RNA tertiary structure: RNA motifs, patterns, and networks. ACC Chem. Res.<\/i> 44<\/b>, 1302\u20131311 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  2. \n

    Klein, D. J., Schmeing, T. M., Moore, P. B. & Steitz, T. A. The kink-turn: a new RNA secondary structure motif. EMBO J.<\/i> 20<\/b>, 4214\u20134221 (2001).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  3. \n

    Lilley, D. M. The K-turn motif in riboswitches and other RNA species. Biochim. Biophys. Acta<\/i> 1839<\/b>, 995\u20131004 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  4. \n

    Schroeder, K. T., McPhee, S. A., Ouellet, J. & Lilley, D. M. A structural database for k-turn motifs in RNA. RNA<\/i> 16<\/b>, 1463\u20131468 (2010).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  5. \n

    Rozhdestvensky, T. S. et al. Binding of L7Ae protein to the K-turn of archaeal snoRNAs: a shared RNA binding motif for C\/D and H\/ACA box snoRNAs in Archaea. Nucleic Acids Res.<\/i> 31<\/b>, 869\u2013877 (2003).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  6. \n

    Baird, N. J., Zhang, J., Hamma, T. & Ferr\u00e9-D\u2019Amar\u00e9, A. R. YbxF and YlxQ are bacterial homologs of L7Ae and bind K-turns but not K-loops. RNA<\/i> 18<\/b>, 759\u2013770 (2012).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  7. \n

    Nottrott, S. et al. Functional interaction of a novel 15.5kD [U4\/U6.U5] tri-snRNP protein with the 5\u2032 stem-loop of U4 snRNA. EMBO J.<\/i> 18<\/b>, 6119\u20136133 (1999).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  8. \n

    Szewczak, L. B. W., DeGregorio, S. J., Strobel, S. A. & Steitz, J. A. Exclusive interaction of the 15.5 kD protein with the terminal box C\/D motif of a methylation guide snoRNP. Chem. Biol.<\/i> 9<\/b>, 1095\u20131107 (2002).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  9. \n

    Szewczak, L. B. W., Gabrielsen, J. S., Degregorio, S. J., Strobel, S. A. & Steitz, J. A. Molecular basis for RNA kink-turn recognition by the h15.5K small RNP protein. RNA<\/i> 11<\/b>, 1407\u20131419 (2005).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  10. \n

    Chawla, M., Oliva, R., Bujnicki, J. M. & Cavallo, L. An atlas of RNA base pairs involving modified nucleobases with optimal geometries and accurate energies. Nucleic Acids Res.<\/i> 43<\/b>, 6714\u20136729 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  11. \n

    Polikanov, Y. S., Melnikov, S. V., S\u00f6ll, D. & Steitz, T. A. Structural insights into the role of rRNA modifications in protein synthesis and ribosome assembly. Nat. Struct. Mol. Biol.<\/i> 22<\/b>, 342\u2013344 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  12. \n

    Liu, S., Ghalei, H., L\u00fchrmann, R. & Wahl, M. C. Structural basis for the dual U4 and U4atac snRNA-binding specificity of spliceosomal protein hPrp31. RNA<\/i> 17<\/b>, 1655\u20131663 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  13. \n

    Vidovic, I., Nottrott, S., Hartmuth, K., L\u00fchrmann, R. & Ficner, R. Crystal structure of the spliceosomal 15.5kD protein bound to a U4 snRNA fragment. Mol. Cell<\/i> 6<\/b>, 1331\u20131342 (2000).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  14. \n

    Jafarifar, F., Dietrich, R. C., Hiznay, J. M. & Padgett, R. A. Biochemical defects in minor spliceosome function in the developmental disorder MOPD I. RNA<\/i> 20<\/b>, 1078\u20131089 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  15. \n

    Edery, P. et al. Association of TALS developmental disorder with defect in minor splicing component U4atac snRNA. Science<\/i> 332<\/b>, 240\u2013243 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  16. \n

    He, H. et al. Mutations in U4atac snRNA, a component of the minor spliceosome, in the developmental disorder MOPD I. Science<\/i> 332<\/b>, 238\u2013240 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  17. \n

    Merico, D. et al. Compound heterozygous mutations in the noncoding RNU4ATAC<\/i> cause Roifman Syndrome by disrupting minor intron splicing. Nat. Commun.<\/i> 6<\/b>, 8718 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  18. \n

    Farach, L. S. et al. The expanding phenotype of RNU4ATAC<\/i> pathogenic variants to Lowry Wood syndrome. Am. J. Med. Genet. A<\/i> 176<\/b>, 465\u2013469 (2018).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  19. \n

    Lapinaite, A. et al. The structure of the box C\/D enzyme reveals regulation of RNA methylation. Nature<\/i> 502<\/b>, 519\u2013523 (2013).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  20. \n

    Siegfried, N. A., Busan, S., Rice, G. M., Nelson, J. A. E. & Weeks, K. M. RNA motif discovery by SHAPE and mutational profiling (SHAPE-MaP). Nat. Methods<\/i> 11<\/b>, 959\u2013965 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  21. \n

    Luo, Q.-J. et al. RNA structure probing reveals the structural basis of Dicer binding and cleavage. Nat. Commun.<\/i> 12<\/b>, 3397 (2021).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  22. \n

    McPhee, S. A., Huang, L. & Lilley, D. M. A critical base pair in k-turns that confers folding characteristics and correlates with biological function. Nat. Commun.<\/i> 5<\/b>, 5127 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  23. \n

    Huang, L. et al. Structure and folding of four putative kink turns identified in structured RNA species in a test of structural prediction rules. Nucleic Acids Res.<\/i> 49<\/b>, 5916\u20135924 (2021).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  24. \n

    Liu, J. & Lilley, D. M. The role of specific 2\u2032-hydroxyl groups in the stabilization of the folded conformation of kink-turn RNA. RNA<\/i> 13<\/b>, 200\u2013210 (2007).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  25. \n

    Huang, L., Wang, J. & Lilley, D. M. A critical base pair in k-turns determines the conformational class adopted, and correlates with biological function. Nucleic Acids Res.<\/i> 44<\/b>, 5390\u20135398 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  26. \n

    Ashraf, S., Huang, L. & Lilley, D. M. J. Effect of methylation of adenine N6<\/sup> on kink turn structure depends on location. RNA Biol.<\/i> 16<\/b>, 1377\u20131385 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  27. \n

    Huang, L., Ashraf, S., Wang, J. & Lilley, D. M. Control of box C\/D snoRNP assembly by N6<\/sup>-methylation of adenine. EMBO Rep.<\/i> 18<\/b>, 1631\u20131645 (2017).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  28. \n

    Xuan, J. J. et al. RMBase v2.0: deciphering the map of RNA modifications from epitranscriptome sequencing data. Nucleic Acids Res.<\/i> 46<\/b>, D327\u2013D334 (2018).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  29. \n

    Zhipeng, L. et al. RNA duplex map in living cells reveals higher-order transcriptome structure. Cell<\/i> 165<\/b>, 1267\u20131279 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  30. \n

    Tycowski, K. T., Aab, A. & Steitz, J. A. Guide RNAs with 5\u2032 caps and novel box C\/D snoRNA-like domains for modification of snRNAs in metazoa. Curr. Biol.<\/i> 14<\/b>, 1985\u20131995 (2004).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  31. \n

    Cologne, A. et al. New insights into minor splicing-a transcriptomic analysis of cells derived from TALS patients. RNA<\/i> 25<\/b>, 1130\u20131149 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  32. \n

    Will, C. L. et al. The human 18S U11\/U12 snRNP contains a set of novel proteins not found in the U2-dependent spliceosome. RNA<\/i> 10<\/b>, 929\u2013941 (2004).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  33. \n

    Dominguez, D. et al. Sequence, structure, and context preferences of human RNA binding proteins. Mol. Cell<\/i> 70<\/b>, 854\u2013867 e859 (2018).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  34. \n

    Jolma, A. et al. Binding specificities of human RNA-binding proteins toward structured and linear RNA sequences. Genome Res.<\/i> 30<\/b>, 962\u2013973 (2020).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  35. \n

    Engreitz, J. M. et al. Local regulation of gene expression by lncRNA promoters, transcription and splicing. Nature<\/i> 539<\/b>, 452\u2013455 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  36. \n

    Anzalone, A. V. et al. Search-and-replace genome editing without double-strand breaks or donor DNA. Nature<\/i> 576<\/b>, 149\u2013157 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  37. \n

    Alioto, T. S. U12DB: a database of orthologous U12-type spliceosomal introns. Nucleic Acids Res.<\/i> 35<\/b>, D110\u2013D115 (2007).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  38. \n

    Olthof, A. M., Hyatt, K. C. & Kanadia, R. N. Minor intron splicing revisited: identification of new minor intron-containing genes and tissue-dependent retention and alternative splicing of minor introns. BMC Genomics<\/i> 20<\/b>, 686 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  39. \n

    Madan, V. et al. Aberrant splicing of U12-type introns is the hallmark of ZRSR2 mutant myelodysplastic syndrome. Nat. Commun.<\/i> 6<\/b>, 6042 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  40. \n

    Reber, S. et al. Minor intron splicing is regulated by FUS and affected by ALS-associated FUS mutants. EMBO J.<\/i> 35<\/b>, 1504\u20131521 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  41. \n

    Verberne, E. A., Faries, S., Mannens, M., Postma, A. V. & van Haelst, M. M. Expanding the phenotype of biallelic RNPC3<\/i> variants associated with growth hormone deficiency. Am. J. Med. Genet. A<\/i> 182<\/b>, 1952\u20131956 (2020).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  42. \n

    Argente, J. et al. Defective minor spliceosome mRNA processing results in isolated familial growth hormone deficiency. EMBO Mol. Med.<\/i> 6<\/b>, 299\u2013306 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  43. \n

    Martos-Moreno, G. et al. Response to growth hormone in patients with RNPC3<\/i> mutations. EMBO Mol. Med.<\/i> 10<\/b>, e9143 (2018).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  44. \n

    Elsaid, M. F. et al. Mutation in noncoding RNA RNU12 causes early onset cerebellar ataxia. Ann. Neurol.<\/i> 81<\/b>, 68\u201378 (2017).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  45. \n

    Burns, R. et al. Homozygous splice mutation in CWF19L1<\/i> in a Turkish family with recessive ataxia syndrome. Neurology<\/i> 83<\/b>, 2175\u20132182 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  46. \n

    Evers, C. et al. Exome sequencing reveals a novel CWF19L1<\/i> mutation associated with intellectual disability and cerebellar atrophy. Am. J. Med. Genet. A<\/i> 170<\/b>, 1502\u20131509 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  47. \n

    Nguyen, M. et al. Pathogenic CWF19L1<\/i> variants as a novel cause of autosomal recessive cerebellar ataxia and atrophy. Eur. J. Hum. Genet.<\/i> 24<\/b>, 619\u2013622 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  48. \n

    Bailey, T. L. & Elkan, C. Fitting a mixture model by expectation maximization to discover motifs in biopolymers. Proc. Int. Conf. Intell. Syst. Mol. Biol.<\/i> 2<\/b>, 28\u201336 (1994).<\/p>\n

    CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  49. \n

    Kramer, K. et al. Photo-cross-linking and high-resolution mass spectrometry for assignment of RNA-binding sites in RNA-binding proteins. Nat. Methods<\/i> 11<\/b>, 1064\u20131070 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  50. \n

    Morlan, J. D., Qu, K. & Sinicropi, D. V. Selective depletion of rRNA enables whole transcriptome profiling of archival fixed tissue. PLoS ONE<\/i> 7<\/b>, e42882 (2012).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  51. \n

    Adiconis, X. et al. Comparative analysis of RNA sequencing methods for degraded or low-input samples. Nat. Methods<\/i> 10<\/b>, 623\u2013629 (2013).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  52. \n

    Zinshteyn, B., Wangen, J. R., Hua, B. & Green, R. Nuclease-mediated depletion biases in ribosome footprint profiling libraries. RNA<\/i> 26<\/b>, 1481\u20131488 (2020).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  53. \n

    Gillen, A. E., Yamamoto, T. M., Kline, E., Hesselberth, J. R. & Kabos, P. Improvements to the HITS-CLIP protocol eliminate widespread mispriming artifacts. BMC Genomics<\/i> 17<\/b>, 338 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  54. \n

    Keene, J. D., Komisarow, J. M. & Friedersdorf, M. B. RIP-Chip: the isolation and identification of mRNAs, microRNAs and protein components of ribonucleoprotein complexes from cell extracts. Nat. Protoc.<\/i> 1<\/b>, 302\u2013307 (2006).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  55. \n

    Deng, B. et al. An LTR retrotransposon-derived lncRNA interacts with RNF169 to promote homologous recombination. EMBO Rep.<\/i> 20<\/b>, e47650 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  56. \n

    Spitale, R. C. et al. Structural imprints in vivo decode RNA regulatory mechanisms. Nature<\/i> 519<\/b>, 486\u2013490 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  57. \n

    Phelps, W. A., Carlson, A. E. & Lee, M. T. Optimized design of antisense oligomers for targeted rRNA depletion. Nucleic Acids Res.<\/i> 49<\/b>, e5 (2021).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  58. \n

    Chu, C., Qu, K., Zhong, F. L., Artandi, S. E. & Chang, H. Y. Genomic maps of long noncoding RNA occupancy reveal principles of RNA-chromatin interactions. Mol. Cell<\/i> 44<\/b>, 667\u2013678 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  59. \n

    Percharde, M. et al. A LINE1-nucleolin partnership regulates early development and ESC identity. Cell<\/i> 174<\/b>, 391\u2013405 (2018).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  60. \n

    Chu, C. et al. Systematic discovery of Xist RNA binding proteins. Cell<\/i> 161<\/b>, 404\u2013416 (2015).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  61. \n

    Haeussler, M. et al. The UCSC Genome Browser database: 2019 update. Nucleic Acids Res.<\/i> 47<\/b>, D853\u2013d858 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  62. \n

    Harrow, J. et al. GENCODE: the reference human genome annotation for The ENCODE Project. Genome Res.<\/i> 22<\/b>, 1760\u20131774 (2012).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  63. \n

    Lestrade, L. & Weber, M. J. snoRNA-LBME-db, a comprehensive database of human H\/ACA and C\/D box snoRNAs. Nucleic Acids Res.<\/i> 34<\/b>, D158\u2013D162 (2006).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  64. \n

    Yang, J. H., Shao, P., Zhou, H., Chen, Y. Q. & Qu, L. H. deepBase: a database for deeply annotating and mining deep sequencing data. Nucleic Acids Res.<\/i> 38<\/b>, D123\u2013D130 (2010).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  65. \n

    Jorjani, H. et al. An updated human snoRNAome. Nucleic Acids Res.<\/i> 44<\/b>, 5068\u20135082 (2016).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  66. \n

    Kishore, S. et al. Insights into snoRNA biogenesis and processing from PAR-CLIP of snoRNA core proteins and small RNA sequencing. Genome Biol.<\/i> 14<\/b>, R45 (2013).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  67. \n

    Pruitt, K. D. et al. RefSeq: an update on mammalian reference sequences. Nucleic Acids Res.<\/i> 42<\/b>, D756\u2013D763 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  68. \n

    Quinlan, A. R. & Hall, I. M. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics<\/i> 26<\/b>, 841\u2013842 (2010).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  69. \n

    Camacho, C. et al. BLAST+: architecture and applications. BMC Bioinformatics<\/i> 10<\/b>, 421 (2009).<\/p>\n

    Article<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  70. \n

    Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J.<\/i> 17<\/b>, 10\u201312 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  71. \n

    Magoc, T. & Salzberg, S. L. FLASH: fast length adjustment of short reads to improve genome assemblies. Bioinformatics<\/i> 27<\/b>, 2957\u20132963 (2011).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  72. \n

    Dobin, A. et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics<\/i> 29<\/b>, 15\u201321 (2013).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  73. \n

    Frankish, A. et al. GENCODE reference annotation for the human and mouse genomes. Nucleic Acids Res.<\/i> 47<\/b>, D766\u2013D773 (2019).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n

  74. \n

    Shen, S. et al. rMATS: robust and flexible detection of differential alternative splicing from replicate RNA-seq data. Proc. Natl Acad. Sci. USA<\/i> 111<\/b>, E5593\u2013E5601 (2014).<\/p>\n

    Article<\/a>\u00a0
    \n     CAS<\/a>\u00a0
    \n     PubMed<\/a>\u00a0
    \n     PubMed Central<\/a>\u00a0
    \n
    \n Google Scholar<\/a>\u00a0\n <\/p>\n<\/li>\n<\/ol>\n

    Download references<\/a><\/p>\n<\/div>\n<\/div>\n

    \n

    Acknowledgements<\/h2>\n

    We thank L. Huang from Sun Yat-sen Memorial Hospital for valuable suggestions on this manuscript. We thank Y. Zhang for sharing Sleeping Beauty transposon system. We thank all the staff from Sun Yat-sen University for their support and hard work during the COVID-19 pandemic. This work was supported, in part, by the National Key R&D Program of China (2019YFA0802202 (to J.Y.) and 2022YFA1303300 (to J.Y.)); the National Natural Science Foundation of China (32225011 (to J.Y.), 91940304 (to J.Y.), 31971228 (to J.Y.), 31770879 (to J.Y.), 31970604 (to L.Q.), 31900903 (to B.L.) and 32100467(to S.L.)); the Youth Science and Technology Innovation Talent of Guangdong TeZhi Plan (2019TQ05Y181 (to J.Y.)); funds from Guangzhou City (202002030351 (to J.Y.)); and Fundamental Research Funds for the Central Universities, Sun Yat-sen University (20lgpy112 (to B.L.) and 2021qntd26 (to B.L.)).<\/p>\n<\/div>\n

    \n

    Author information<\/h2>\n
    \n

    Author notes<\/span><\/p>\n

      \n
    1. \n

      These authors contributed equally: Bin Li, Shurong Liu, Wujian Zheng, Anrui Liu, Peng Yu.<\/p>\n<\/li>\n<\/ol>\n

      Authors and Affiliations<\/h3>\n
        \n
      1. \n

        MOE Key Laboratory of Gene Function and Regulation, State Key Laboratory of Biocontrol, School of Life Sciences, Sun Yat-sen University, Guangzhou, China<\/p>\n

        Bin Li,\u00a0Shurong Liu,\u00a0Wujian Zheng,\u00a0Anrui Liu,\u00a0Peng Yu,\u00a0Di Wu,\u00a0Jie Zhou,\u00a0Ping Zhang,\u00a0Chang Liu,\u00a0Qiao Lin,\u00a0Jiayi Ye,\u00a0Qiaojuan Huang,\u00a0Hui Zhou,\u00a0Lianghu Qu\u00a0&\u00a0Jianhua Yang<\/p>\n<\/li>\n

      2. \n

        The Fifth Affiliated Hospital, Sun Yat-sen University, Zhuhai, China<\/p>\n

        Wujian Zheng,\u00a0Anrui Liu,\u00a0Di Wu,\u00a0Jie Zhou,\u00a0Qiao Lin,\u00a0Jiayi Ye,\u00a0Simeng He\u00a0&\u00a0Jianhua Yang<\/p>\n<\/li>\n

      3. \n

        Department of Systems Biology, Beckman Research Institute of City of Hope, Monrovia, CA, USA<\/p>\n

        Jianjun Chen<\/p>\n<\/li>\n<\/ol>\n

        Contributions<\/h3>\n

        J.Y., B.L. and L.Q. conceived and designed the entire project. J.Y. and L.Q. designed and supervised the research. B.L., S.L., W.Z., A.L., P.Y., D.W., J.Z., P.Z., C.L., Q.L., J.Y., S.H., Q.H., H.Z. and J.Y. performed the experiments and\/or data analyses. J.Y., B.L. and A.L. performed the genome-wide or transcriptome-wide data analyses. J.Y., B.L., S.L. and L.Q. contributed reagents\/analytic tools and\/or grant support. J.C. provided helpful discussions. J.Y., L.Q., B.L., S.L. and A.L. wrote and revised the paper. All authors discussed the results and commented on the manuscript.<\/p>\n

        Corresponding authors<\/h3>\n

        Correspondence to
        \n Lianghu Qu<\/a> or Jianhua Yang<\/a>.<\/p>\n<\/div>\n<\/div>\n

        \n

        Ethics declarations<\/h2>\n
        \n

        Competing interests<\/h3>\n

        J.C. is a scientific advisory board member of Race Oncology. The remaining authors declare no competing interests.<\/p>\n<\/p><\/div>\n<\/div>\n

        \n

        Peer review<\/h2>\n
        \n

        Peer review information<\/h3>\n

        Nature Biotechnology<\/i> thanks the anonymous reviewers for their contribution to the peer review of this work.<\/p>\n<\/p><\/div>\n<\/div>\n

        \n

        Additional information<\/h2>\n

        Publisher\u2019s note<\/b> Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.<\/p>\n<\/div>\n

        \n

        Extended data<\/h2>\n
        \n
        \n

        Extended Data Fig. 1 Characterization of the RIP-PEN-seq technique.<\/a><\/h3>\n
        \n

        a<\/b>, Diagram of RNase H-based rRNA depletion for the construction of RIP-PEN-seq library. b<\/b>, Western blotting analysis of the overexpression of 15.5K-FLAG protein in HEK293T cells. GAPDH serves as the loading control. c<\/b>, RIP from HEK293T cells was performed using anti-FLAG and IgG. d<\/b>, e<\/b>, Meta-analyses of the RIP-PEN-seq results for the start (d<\/b>) and end (e<\/b>) sites of the forward K-turn RNAs (box C\/D ncRNAs) in HEK293T cells. f<\/b>, Genome-browser plot of RIP-PEN-seq (coverage, green; 5\u2019-start, red; 3\u2019-end, purple) for representative forward K-turn RNAs (box C\/D snoRNAs) in the introns of GAS5. g<\/b>, Computational workflow for analysis of the RIP-PEN-seq sequencing data and identification of candidate transcripts. h<\/b>, KturnSeeker core algorithm workflow. KturnSeeker was developed to identify and quantify forward (fktRNAs) and backward ktRNAs (bktRNAs) from RIP-PEN-seq data. KturnSeeker can screen bktRNAs as well as fktRNAs by reverse searching the K-turn structure. i<\/b>, Gene type distribution of forward ktRNAs identified from RIP-PEN-seq. CD-type represents C\/D box-containing snoRNAs and scaRNAs. ACA-type represents snoRNAs or scaRNAs that only contain H\/ACA boxes.<\/p>\n

        \n Source data<\/a>\n <\/p>\n<\/div>\n<\/div>\n

        \n

        Extended Data Fig. 2 Structural characterization of bktRNAs.<\/a><\/h3>\n

        a<\/b>, The position and significance of motifs located within bktRNAs. The enriched motifs were identified by MEME software. The K-turn structural motif of bktRNAs is composed of two conserved sequence motifs: a CUGA motif often 4 nt downstream of the 5\u2019 end and a UGAUG motif often 2 nt upstream of the 3\u2019 end. The position p-value is defined as the probability that a random sequence would have a motif match score greater or equal to the sequence under test. b<\/b>, Another twelve novel bktRNAs discovered from RIP-PEN-seq data. All novel bktRNAs had a CUGA motif that was often 4 nt downstream of the 5\u2019 end and a UGAUG motif that was often 2 nt upstream of the 3\u2019 end. These two sequence motifs were located within the K-turn structural motifs of bktRNAs. The 5\u2019 motif (CUGA) and 3\u2019 motif (UGAUG) are marked with black rectangles. The non-canonical A\u2022G and G\u2022A base pairs and the mismatch in the backward K-turn structure are also marked. c<\/b>, Schematic overview of 15.5\u2009K RIP-PEN-SHAPE-MaP. d<\/b>, Secondary structure of consensus fktRNA (left panel, n\u2009=\u200998). Violin plots displaying the SHAPE reactivity across the forward K-turn structure (including 5\u2019 C stem, internal (internal loop), 5\u2019 NC stem, Loop, 3\u2019 NC stem, and 3\u2019 C stem), averaged across all known box C\/D snoRNAs (right panel). The boxplots indicate the median and the upper and lower quartiles. e<\/b>, The predicted secondary structure (upper panel) and SHAPE reactivity signal (lower panel) on fktRNA SNORD102 (also known as U102). The forward K-turn structure is indicated in the structure figure. The NC-stem and C-stem are marked with black and red underlines in the bar plot, respectively. The SHAPE reactivity signal was determined by RIP-PEN-SHAPE-MaP in this study.<\/p>\n<\/div>\n

        \n

        Extended Data Fig. 3 The backward K-turn sequence composition, predicted functions and tissue-specific expression profile of bktRNAs.<\/a><\/h3>\n

        a<\/b>, The secondary structure of consensus forward K-turn RNA (fktRNA) and backward K-turn RNA (bktRNA). The nucleotide positions in the K-turn structure are named according to the nomenclature rules for the forward K-turn structure. b<\/b>, Matrix plot showing the number of human bktRNAs with the indicated nucleotide in the 3b:3n sequences. c<\/b>, Number of human bktRNAs with the four possible Watson-Crick base pairs in the -1b:-1n position. d<\/b>, Matrix plot showing the number of mouse bktRNAs with the indicated nucleotide in the 3b:3n sequences. e<\/b>, Number of mouse bktRNAs with the four possible Watson-Crick base pairs in the -1b:-1n position. f<\/b>, Number of bktRNAs with or without m6<\/sup>A modification in humans and mice. g<\/b>, Number of bktRNAs with or without m6<\/sup>A modification at the 1n position in humans and mice. h<\/b>, Enrichment analysis of the bktRNA host protein-coding genes by Metascape software. i<\/b>, Tissue-specific expression profiles of bktRNAs. The expression levels of bktRNAs are displayed in the rows and the tissues are shown in the columns. The rows and columns are sorted based on k-means clustering of bktRNA genes. The colour intensity represents the tissue-specific score (JS score) as calculated for each bktRNA using the csSpecificity function. Representative bktRNAs are indicated in the right panel.<\/p>\n<\/div>\n

        \n

        Extended Data Fig. 4 Genomic characterization, expression, conservation, and secondary structure of bktRNA1.<\/a><\/h3>\n

        a<\/b>, Genome-browser plot of RIP-PEN-seq (coverage, blue; 5\u2019-start, red; 3\u2019-end, yellow) for bktRNA1, as well as the evolutionary conservation across 100 vertebrates (green). b<\/b>, Secondary structure of bktRNA1 in the human genome was predicted by R-scape software. The SHAPE reactivities for each nucleotide were mapped to secondary structures using R2R software. The box H\/ACA domain is indicated with a black dashed box, and the backward K-turn structure and the potential K-turn-like structure are marked with green dashed boxes. The NC stem and C stem are indicated with black lines. The blue boxes show the representative motifs. CAB, Cajal body box. c<\/b>, The SHAPE reactivity signal on bktRNA1. The representative motifs are underlined in the bar plot. The SHAPE reactivity signal was determined by RIP-PEN-SHAPE-MaP in this study.<\/p>\n<\/div>\n

        \n

        Extended Data Fig. 5 Secondary structures and subcellular localization of bktRNA1 and its interacting partner U12 snRNA.<\/a><\/h3>\n

        a<\/b>, Predicted conserved RNA structure of bktRNA1 determined by measuring pairwise covariations with R-scape software. The H\/ACA domain is indicated with a black dashed box. The functional region paired with U12 snRNA is indicated with a blue dashed box. b<\/b>, In situ co-localization of bktRNA1 with 15.5\u2009K proteins and U12 snRNAs in HEK293T cells by fluorescent in situ hybridization (FISH) and immunofluorescence (IF) microscopy. White arrows indicate the signal detected by probes or antibody. c<\/b>, In situ co-localization of bktRNA1 with 15.5\u2009K proteins and U12 snRNAs in HCT116 cells by fluorescent in situ hybridization (FISH) and immunofluorescence (IF) microscopy. White arrows indicate the signal detected by probes or antibody.<\/p>\n<\/div>\n

        \n

        Extended Data Fig. 6 Splicing efficiency analysis for wild-type and knockout bktRNA1.<\/a><\/h3>\n

        a<\/b>, Workflow for intron retention analysis in HCT116 and KO-bktRNA1 cells. b<\/b>, Proportion of aberrantly retained U12- and U2-type introns (filtered by p\u2009<\u20090.05) in bktRNA1-deficient cells. c<\/b>, Proportion of statistically significant changes (filtered by p\u2009<\u20090.05) in U12- and U2-type genes in bktRNA1-deficient cells. d<\/b>, The ratio of spliced to unspliced pre-mRNA for U12-type introns was determined by qPCR in bktRNA1-deficient cells. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. e<\/b>, The ratio of spliced to unspliced pre-mRNA for U12-type and U2-type (GAPDH) introns was determined by qPCR in bktRNA1-rescued HCT116 KO-4 cells. f<\/b>, The ratio of spliced to unspliced pre-mRNA for U12-type and U2-type (GAPDH) introns was determined by qPCR in SNORA12-rescued HCT116 KO-4 cells. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. ns, no significance. g<\/b>, The ratio of spliced to unspliced pre-mRNA for U12-type and U2-type (GAPDH) introns was determined by qPCR in artificial scaRNA-overexpressing HCT116 KO-4 cells. Data are presented as mean values +\/- SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. ns, no significance.<\/p>\n<\/div>\n

        \n

        Extended Data Fig. 7 Depletion of bktRNA1 affects the interaction between U12 and ZCRB1.<\/a><\/h3>\n
        \n

        a<\/b>, Western blots showing precipitation with each indicated antibody in wild-type (WT) and bktRNA1-deficient KO-4 (KO) cells. b<\/b>, Native RIP was performed in wild-type (WT) and bktRNA1-deficient KO-4 (KO) cells using each indicated antibody or normal IgG antibody, after which qPCR was performed with primers recognizing minor splice snRNAs (U11, U12, U4atac, U5, U6atac). The percentage of RIP-enriched snRNAs was calculated relative to the input RNA. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. ns, no significance. c<\/b>, ZCRB1 RIP-enriched snRNAs were detected by Northern blotting in wild-type (WT) and bktRNA1-deficient KO-4 (KO) cells. U6 snRNA served as a negative control. d<\/b>, ZCRB1 RIP-enriched snRNAs were detected by Northern blotting in bktRNA1-deficient and bktRNA1-rescued cells. U6 snRNA served as a negative control.<\/p>\n

        \n Source data<\/a>\n <\/p>\n<\/div>\n<\/div>\n

        \n

        Extended Data Fig. 8 ZCRB1 knockdown affects U12-type intron splicing.<\/a><\/h3>\n
        \n

        a<\/b>, qPCR (upper panel) and western blotting analysis (lower panel) of Dox-inducible ZCRB1 knockdown in HCT116 cells. GAPDH was used as an internal reference gene for qPCR, and GAPDH served as the loading control for western blotting. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. b<\/b>, c<\/b>, Dot plots displaying the intron retention levels in a representative pairwise analysis of ZCRB1 knockdown and negative control cells. The red dots represent U12-type introns, and the blue dots represent U2-type introns. d<\/b>, Proportion of aberrantly retained U2- and U12- type introns in ZCRB1 knockdown cells. The red boxes represent retained introns, and the blue boxes represent unretained introns. e<\/b>, f<\/b>, Cumulative fraction of the inclusion level difference between U12-type and U2-type introns in ZCRB1 knockdown and negative control cells. The P<\/i> value on the cumulative plots of inclusion level differences were calculated using a two-sided Mann-Whitney-Wilcoxon test. g<\/b>, Venn diagram showing the numbers of overlapping retained introns across four bktRNA1-deficient HCT116 cell lines and ZCRB1 knockdown cells. h<\/b>, CCK-8 assay of HCT116 cells with bktRNA1 knockout. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. i<\/b>, Colony formation assay of HCT116 cells with bktRNA1 knockout. j<\/b>, Quantitative analysis of colony formation assay in the indicated lines. Data are presented as mean values +\/- SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. k<\/b>, CCK-8 assay of HCT116 cells with ZCRB1 knockdown. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test. l<\/b>, Colony formation assay of HCT116 cells with ZCRB1 knockdown. m<\/b>, Quantitative analysis of colony formation assay in the indicated lines. Data are presented as mean values +\/\u2212 SEM (n\u2009=\u20093, biological replicates), two-tailed, paired t-test.<\/p>\n

        \n Source data<\/a>\n <\/p>\n<\/div>\n<\/div>\n<\/div>\n<\/div>\n

        \n

        Supplementary information<\/h2>\n<\/div>\n
        \n

        Source data<\/h2>\n<\/div>\n