Quantitative analysis of small RNA pseudouridylation reveals interplay of PUS enzymes in tRNA anticodon stem-loop

Domain: neuroscience
Published: 2026-02-16
Editorial outcome: major_revision

Abstract
Pseudouridine (Ψ) is an abundant modification in small RNA catalyzed by multiple pseudouridine synthases (PUSs). However, the substrate specificity of human PUSs remains elusive. Here, we adopted PRAISE, a quantitative Ψ detection method, to profile pseudouridylation in small RNA, including cytosolic and mitochondrial tRNAs, snRNA, and snoRNA. We found that snoRNA pseudouridylation is mediated not only by RNA-guided DKC1, but also by the stand-alone enzyme PUS7 at a specific site. Interestingly, several PUS enzymes, including PUS1, RPUSD1, and PUS7, which install nearby Ψ sites within tRNA anticodon stem-loop, can influence pseudouridylation catalyzed by other PUSs, revealing an unrecognized interplay during Ψ formation. For the three RluA family enzymes, RPUSD1 catalyzes the canonical Ψ30 in tRNA-Ile and Ψ72 in tRNA-Arg isoacceptors. RPUSD2 pseudouridylates Ψ31 of mt-tRNALeu(CUN), Ψ32 of mt-tRNAPro and mt-tRNACys, whereas RPUSD3 lacks tRNA activity. Together, our quantitative Ψ profiling characterized PUS tRNA substrates and revealed unexpected PUS interplay.

Introduction
Pseudouridine (Ψ), known as the “fifth nucleotide” of RNA, is the most abundant post-transcriptional modification and is prevalent in various RNA species1,2. Ψ has unique physiochemical properties due to the additional hydrogen bond donor at the N1 position, which has been shown to increase base-pairing stability3. Emerging evidence suggests that Ψ is dynamically regulated in response to cellular stress and is involved in multiple processes that influence pre-mRNA splicing and translational regulation4,5. Ψ is highly enriched in noncoding RNAs, particularly small RNAs such as tRNA, snRNA, and snoRNA, and is clustered in evolutionarily conserved and functionally important regions6. More than 14 positions have been identified in cytosolic tRNA, which are distributed across critical regions including the D loop, anticodon stem-loop, TΨC loop, and acceptor stem2,7. The widespread occurrence of Ψ in tRNA highlights its regulatory role in tRNA processing and translational control4. Ψ installation is mediated by pseudouridine synthases (PUSs). In humans, 13 PUS enzymes can be classified into five families: RluA, TruA, TruB, TruD, and PUS108. Mutations and dysregulated expression of several PUS enzymes lead to aberrant pseudouridylation of tRNAs, which affects various cellular pathways and contributes to disease progression5. Genetic mutations in PUS1, which is responsible for the pseudouridylation in specific tRNAs9, lead to mitochondrial myopathy and sideroblastic anemia (MLASA)10,11. A mutation in PUS1 in a patient with an MLASA phenotype impairs mitochondrial function and blocks erythroid differentiation12. PUS7 variants are linked to intellectual disability13, and PUS7-dependent Ψ8 modification in tRNA fragments was found to impact translation and contribute to tumorigenesis14. High expression of PUS7 in glioblastoma and its regulation of select tRNA Ψ50 modification are critical for codon-specific translational control of key regulators of GSCs15. Although emerging evidence demonstrates the regulatory roles of human PUSs in physiological and pathological conditions, the tRNA substrates of most human PUS enzymes remain poorly understood. Further investigation of additional PUS targets is crucial for elucidating the regulatory mechanisms of PUSs.

The distribution of Ψ modifications in RNA substrates relies on Ψ detection technologies. Early transcriptome-wide profiling approaches utilized N-Cyclohexyl-N′-(2-morpholinoethyl)carbodiimide methyl-p-toluenesulfonate (CMC) or azide-modified CMC derivatives to label RNA, leveraging CMC-induced reverse transcription termination to identify Ψ sites16–20. DM-Ψ-seq adopts AlkB demethylases to remove methylations that block reverse transcription, thereby improving the sensitivity for Ψ detection in small RNA21. However, it cannot provide quantitative stoichiometric information due to partial labeling efficiency. An nanopore-based technique, Nano-tRNA seq, reported quantitative estimates of tRNA Ψ modification using base-calling errors, and it has been applied to Saccharomyces cerevisiae tRNAs22. We and others have developed quantitative methods such as PRAISE, BID-seq, and BACS, which employ bisulfite and cyclization chemistries to enable quantification of Ψ stoichiometry in the human transcriptome23–25. These quantitative methods facilitate the investigation of the dense Ψ landscape of small RNA substrates mediated by multiple human PUSs.

In this study, we adopted PRAISE, a quantitative Ψ detection method developed by our group, to generate a global map of Ψ modifications in small RNAs of HEK293T cells. We found that pseudouridylation of snoRNAs is mediated not only by RNA-guided enzyme DKC1 but also by the stand-alone enzyme PUS7. By perturbing ten PUS enzymes using stable knockdown or knockout cell lines, including DKC1, PUS1, PUS7, TRUB1, RPUSD1, RPUSD2, RPUSD3, PUS1L, PUS7L, and PUS10, we identified their individual small RNA substrates. Notably, we report for the first time the tRNA targets of the human RluA family and found that the substrate of RPUSD1 differs from that of its yeast counterparts. In addition, our data reveal that PUS enzymes targeting Ψ sites can mutually regulate one another within the tRNA anticodon stem-loop. This includes the interplay between RPUSD1 and PUS1 in the anticodon stem, as well as the mutual influences between PUS7 and PUS1, and between PUS7 and PUS3, in the anticodon loop. Together, we employed quantitative Ψ profiling to identify the small RNA substrates of PUSs and reveal crosstalk among them during tRNA pseudouridylation.

Results
Adopting a quantitative detection method for profiling Ψ sites in small RNA
To achieve quantitative Ψ detection in small RNA less than 200 nt, we sought to implement the optimized bisulfite chemistry according to our recently published technology, named PRAISE23. However, specific modifications, including m1A and m3C modifications, can block reverse transcription. Inspired by DM-tRNA-seq and DM-Ψ-seq, we aimed to remove these modifications in vitro to improve cDNA synthesis efficiency21,26. We isolated and treated small RNA with AlkB demethylase, and demonstrated that m1A modification nearly disappeared after AlkB treatment using LC-MS detection (Supplementary Fig. 1a). Further integrating ligation-based strand-specific library construction, we established a methodological pipeline for Ψ detection in small RNA (Supplementary Fig. 1b).

We first analyzed read coverage across small RNA transcripts and did not find bias for the 5' or 3' end in the small RNA library (Supplementary Fig. 1c, d). The fractions of different RNA types showed that cytosolic tRNAs (cyto-tRNAs) accounted for 82.15%, mitochondrial tRNAs (mt-tRNAs) for 0.4%, 5S rRNA for 12.2%, and snRNA and snoRNA together for 5.0%, respectively (Fig. 1a). We then employed PRAISE to identify Ψ sites in small RNAs from HEK293T cells. Since cytosolic tRNA was well-known to contain abundant Ψ sites, we globally identified and quantified Ψ sites in the cyto-tRNAs to test the performance of PRAISE in small RNA. In the cyto-tRNAs, 713 and 714 Ψ sites were detected in two technical replicates, respectively, with 705 overlapping Ψ sites, thus demonstrating high reproducibility of PRAISE (Fig. 1b and Supplementary Data 1). Among the 17 non-overlapping sites, the majority of them (12/17) exhibited low modification levels in another replicate (3–5%) (Supplementary Fig. 1e). We further quantified the Ψ level for the individual site based on the deletion signal, and found a high correlation of Ψ modification levels between two technical replicates (R2=0.99) (Fig. 1c). In addition, Ψ sites in mt-tRNAs can also be detected by PRAISE. We identified 31 and 35 Ψ sites in two technical replicates, respectively, with 31 overlapping sites (Supplementary Data 2). Notably, these four non-overlapping sites that had been previously validated by mass spectrometry data27 also exhibited detectable modification signals in another PRAISE technical replicate (Supplementary Data 2). These results demonstrated the high reproducibility of PRAISE for Ψ detection in small RNA.

We compared PRAISE with DM-Ψ-seq, an existing CMC-based technique for detecting Ψ sites in small RNAs21. In general, 577 Ψ sites in cyto-tRNAs identified by PRAISE can overlap with Ψ sites identified by DM-Ψ-seq (Supplementary Fig. 1f). In particular, PRAISE showed its ability to detect Ψ sites in the 3' end of tRNA transcripts due to its ability to enable reverse transcription read through. As shown in Arg-TCG-5 tRNA, DM-Ψ-seq failed to identify Ψ72, whereas PRAISE successfully detected it and revealed a high modification level of 90.3% (Fig. 1d). Additionally, the modification levels of the 577 overlapping sites detected by DM-Ψ-seq and PRAISE showed a low correlation, which was attributed to the lack of quantitative capacity of DM-Ψ-seq (Supplementary Fig. 1g). We also analyzed the total Ψ levels in cyto-tRNA detected by PRAISE and BACS and found that both methods demonstrated good performance in quantitative Ψ detection, comparable to the Ψ/U ratio measured by quantitative LC-MS/MS (Supplementary Fig. 1h). We observed a high overlap ratio, with 510 of 542 (94.1%) BACS sites detected by PRAISE (Supplementary Fig. 1i). For 32 BACS-only Ψ sites, the majority of them exhibited deletion rates greater than 2% but below our 5% cutoff (Supplementary Fig. 1j). In addition, the Ψ levels quantified by PRAISE in HEK293T cells are correlated with those reported by BACS in HeLa cells (Supplementary Fig. 1k). These data indicated that PRAISE is suitable for Ψ detection and quantification in small RNA.

Quantitative Ψ landscapes of cytosolic and mitochondrial tRNA
We next aimed to provide a comprehensive quantitative landscape of Ψ modification in both cyto-tRNAs and mt-tRNAs in HEK293T cells. We first analyzed the distribution and levels of 705 Ψ sites in cyto-tRNAs. Using PRAISE for quantitative profiling at the isodecoder level, we found distinct Ψ distribution patterns across tRNA isodecoders (Fig. 1e and Supplementary Data 1). Although Ψ is known to be the most common modification in total cytosolic tRNAs, the number of Ψ sites per modified tRNA molecule ranges from one to six (Supplementary Data 3). Human Leu-AAG-4 tRNA has the highest number of Ψ sites (six sites). When analyzing Ψ sites in different tRNA regions, we found that Ψ occurs more frequently in the anticodon stem loop (332 sites) compared to the TΨC loop (250 sites) and D loop (90 sites), suggesting that Ψ plays a critical role in translational decoding. Among all Ψ sites, Ψ55 in the TΨC loop is the most conserved modification, being present in all cyto-tRNAs. Following this, Ψ38-40 and Ψ27-28 in the anticodon stem are present in 148 of 198 and 132 of 202 cyto-tRNAs that contain uridine at these positions, respectively. Ψ13 in the D loop is relatively less conserved and is present in several types of cyto-tRNAs. Additionally, specific Ψ modification sites are found in individual tRNA; for instance, Ψ72 is present only in the acceptor stem of tRNA-Arg, Ψ8 is present only in tRNA-Gly, and Ψ20 is present only in the D loop of tRNA-Leu, indicating their particular functional and physiological roles in tRNA.

We further quantified the modification levels of individual Ψ sites using PRAISE. We integrated and displayed the median modification level for each position across modified cyto-tRNAs, based on an alignment of tRNA sequences28 (Fig. 1f). From the results, we observed that the universally conserved modification site, Ψ55, exhibited a very high modification level, reaching a median of 83.4%. Ψ35, located at the center of tRNA-Tyr anticodon, showed a high modification level of 86.9%, which may be crucial for the recognition of UAG stop codon29. Other highly-modified sites included Ψ72, Ψ13, and Ψ39, with median modification levels of 62.9, 62.1, and 59.8%, respectively. In contrast, Ψ30-32 showed moderate modification levels of 16.7, 24.3, and 24.6%, respectively. Altogether, we provided a global overview of Ψ profiles and stoichiometries in human cyto-tRNAs.

We next analyzed the Ψ distribution and levels in mt-tRNAs and found Ψ sites in 17 out of the total 22 mt-tRNAs (Supplementary Fig. 1l, m). Compared to cyto-tRNAs, mt-tRNAs display a distinct Ψ landscape. Notably, the highly-modified Ψ13 and Ψ72 sites common in cyto-tRNAs are absent, while Ψ68, though present in the acceptor stem, shows a low modification level (7.3%). Instead, Ψ27 and Ψ28 are the most conserved and highly-modified sites in mt-tRNAs, with median modification levels of 79.9 and 68.3%, respectively (Fig. 1f and Supplementary Data 2). Ψ39 is also conserved, but its modification level varies widely—for example, 12.9% in mt-tRNA Phe versus 84.3% in mt-tRNA His. Surprisingly, Ψ55, despite being the most conserved site in cyto-tRNAs, is found in only six mt-tRNAs, though its modification levels in mt-tRNA Leu(UUR) and mt-tRNA Ser are as high as 94.5 and 89.2%, respectively. Thus, these results suggested the presence of a differential Ψ installation machinery in mt-tRNAs.

PUS-dependent and dynamic Ψ sites in snRNA and snoRNA
We then profiled the quantitative Ψ landscape in snRNA and snoRNA from HEK293T cells. Notably, Ψ modification levels exhibit distinct patterns between snRNA and snoRNA. In snRNA, 28 Ψ sites were identified, with modification levels clustering primarily around 15 and 60% (Fig. 2a, b and Supplementary Data 4). These sites are predominantly located near the 5' end of snRNA sequences30 (Fig. 2c). In snoRNAs, we identified 37, 31, and 6 Ψ sites in box H/ACA snoRNAs, box C/D snoRNAs, and scaRNAs, respectively, with modification levels mainly clustering around 15% (Fig. 2a, d and Supplementary Data 5). Interestingly, Ψ sites are enriched within functional pockets where snoRNAs engage their RNA targets, consistent with previous findings25. In box H/ACA snoRNAs, Ψ sites are mainly located within the two hairpin regions, while in box C/D snoRNAs, they are primarily found between the D' and C boxes, and between the C' and D boxes (Fig. 2e, f). The results indicated that Ψ modifications may facilitate snoRNA-guided 2'-O-methylation and pseudouridylation.

DKC1 was reported to install Ψ on rRNA31,32, mRNA23,24, snRNA33, and snoRNA19. We then knocked down DKC1 and performed PRAISE to quantify DKC1 targets in snRNA and snoRNA. We set a threshold of 30.0% decrease in modification levels and identified 12 high-confidence DKC1-dependent sites among 28 Ψ sites in snRNA (Fig. 2g, h and Supplementary Data 4). In snoRNA, three types of snoRNAs showed differential patterns in expression and Ψ levels after DKC1 knockdown. The expression levels of box H/ACA snoRNA and scaRNA decreased significantly, while the expression of box C/D snoRNA remained unchanged (Fig. 2i). This is consistent with a previous study demonstrating that DKC1 deficiency impairs H/ACA snoRNP formation, thereby compromising the stability of box H/ACA snoRNAs19. For scaRNAs, an earlier work has shown that domain-specific mutations in DKC1 can disrupt their 3' end processing34. Notably, the reads coverage in the DKC1 sample was still sufficient for the identification of Ψ sites in snoRNAs, and we identified 22 and 26 DKC1-dependent Ψ sites in box H/ACA and box C/D snoRNAs, respectively (Fig. 2j, k and Supplementary Data 5). The Ψ levels of both box C/D snoRNA and box H/ACA snoRNA decreased significantly, whereas the Ψ levels of scaRNA remained unchanged (Fig. 2j, k).

Only 48 out of 74 Ψ sites in snoRNA are DKC1-dependent (Fig. 2j, k), prompting us to investigate whether stand-alone PUSs also modify snoRNA. In the formation of Ψ sites on mRNA, stand-alone PUS enzymes such as PUS1, PUS7, and TRUB1 serve as the primary modifying enzymes, mediating the majority of Ψ modifications. We speculated that they might have more diverse RNA targets, including snRNA and snoRNAs. To systematically evaluate this possibility, we then performed PRAISE on small RNAs using a comprehensive panel of PUS-deficient HEK293T cell lines, including PUS1−/−, PUS7−/−, TRUB1−/−, RPUSD1−/−, RPUSD2−/−, RPUSD3−/−, PUS1L−/−, PUS7L−/−, and PUS10−/− (details in a later section). In general, DKC1-independent Ψ sites in snRNA and snoRNA remained unaffected in stand-alone PUS knockout cells (Supplementary Fig. 2a, b). But interestingly, we found that one specific site, Ψ47 in SNORD65, was modified by PUS7, corresponding to the canonical “UGUAG” motif (Fig. 2l). Re-analysis of the BACS dataset further confirmed that Ψ47 in SNORD65 is PUS7-dependent (Supplementary Fig. 2c). These results demonstrated the presence of both RNA-dependent and RNA-independent mechanisms for Ψ formation in snoRNA.

To further investigate these DKC1-independent sites, we assessed potential snoRNA-guided Ψ installation and found that most snRNA Ψ sites (12/16) match annotated snoRNA guides in established snoRNA–RNA target databases35,36. Moreover, in-silico predictions based on the folding energies between snoRNAs and target sites19 further revealed that many Ψ sites (22/25) within snoRNAs themselves can be paired with specific H/ACA snoRNAs (Supplementary Fig. 2d and Supplementary Data 6). Together, these results suggested that, except for the single PUS7-dependent site, Ψ installation in snRNAs and snoRNAs is largely mediated by snoRNA-guided modification pathways.

Previous studies revealed stress-induced pseudouridylation of U2 snRNA in yeast37. To assess whether dynamics occur in human cells, we subjected HEK293T cells to heat shock and oxidative stress (H₂O₂ treatment). In snRNA, we observed that Ψ89 in U2 showed a significant increase under heat shock, while other snRNA Ψ sites remained unchanged (Supplementary Fig. 2e). Under H₂O₂ treatment, none of the snRNA Ψ sites exhibited significant changes (Supplementary Fig. 2e). In contrast, snoRNAs displayed more pronounced stress-responsive dynamics: four sites significantly decreased following heat shock and two sites decreased under H₂O₂ treatment (Supplementary Fig. 2f). Notably, Ψ95 in SNORA61 was consistently reduced under both conditions, suggesting that this site may represent a conserved stress-responsive Ψ site in human snoRNAs.

tRNA substrates of PUS1
We further investigated the PUS1 substrates in tRNA based on the PRAISE data of WT and PUS1−/− cell lines. In cyto-tRNA, the modification levels of 132 Ψ sites were significantly reduced or completely lost after PUS1 knockout (Fig. 3a, b and Supplementary Data 1). The contexts of these PUS1-dependent Ψ sites exhibited a weak “USR (S=G or C, R=G or A)” motif (Supplementary Fig. 3a). Consistent with previous findings, Ψ27 and Ψ28 were confirmed as major PUS1-dependent sites in cyto-tRNAs, with their median modification levels dropping from 36.8% in WT to 2.5% in PUS1-/- cell lines (Supplementary Data 1). In mt-tRNA, 14 Ψ modification sites were completely absent after PUS1 knockout (Fig. 3c, d and Supplementary Data 2). We found that PUS1 is responsible for the formation of the Ψ27, Ψ28 and Ψ68 sites in mt-tRNAs, with the median level of Ψ27 and Ψ28 sites decreasing from 70.1% in WT to 0.5% in PUS1−/− and the median level of Ψ68 sites decreasing from 7.3% in WT to 0.7% in PUS1−/− (Supplementary Data 1). Furthermore, our results demonstrate that PUS1 knockout does not affect global tRNA expression (Supplementary Fig. 3b, c). These results demonstrated that human PUS1 is a multifunctional enzyme that targets both cyto- and mt-tRNAs, in accordance with previous studies9,10.

Interplay of PUS1 with additional Ψ synthases
During the analysis of PUS1 substrates, we unexpectedly found that the quantitative results showed an interesting pattern. When the modification levels of Ψ27/Ψ28 in cytosolic tRNAs decreased upon PUS1 knockout, the modification levels of nearby Ψ sites in the anticodon stem-loop significantly increased (Fig. 3e). Specifically, the modification level of Ψ28 on Ile-AAT tRNAs decreased from 47.2 to 3.6%, while the level of Ψ30 on the same tRNAs dramatically increased from 19.3 to 56.1% (Fig. 3f and Supplementary Data 1). Similarly, the modification level of Ψ27 on Met-CAT tRNAs significantly decreased from 43.8 to 7.0%, accompanied by an increase in Ψ31 levels from 24.3 to 48.6%. Notably, in Arg-TCT tRNAs, the Ψ28 level decreased from 53.8 to 11.1%, whereas the Ψ36 level increased markedly from 10.7 to 51.4% upon PUS1 knockout (Supplementary Data 1). In addition, after PUS1 knockout, the modification levels of Ψ39 sites also increased when the modification levels of Ψ27/Ψ28 sites decreased within the same tRNA molecule (Fig. 3e). To exclude artificial interplay caused by inefficient labeling of adjacent Ψ sites, we synthesized a spike-in RNA oligo containing two Ψ sites within an Ile-AAT anticodon stem-loop (ASL)–like sequence (Supplementary Fig. 3d). The expected and measured Ψ levels showed excellent agreement, consistent with the single-Ψ spike-in control (Supplementary Fig. 3e). In addition, BACS data further showed the elevated Ψ30 level in Ile-AAT upon PUS1 knockout, providing orthogonal evidence for the interplay (Supplementary Fig. 3f). Together, these results suggested a common interplay within the anticodon stem-loop exists between PUS1 and additional PUSs.

Given that perturbation of one RNA modification enzyme can influence the expression and activity of others38,39, we first analyzed whether PUS1 knockout affects the gene expression of other PUS enzymes. Among the remaining 12 PUSs, the expression levels of most remained unchanged, and only RPUSD1 and RPUSD3, both members of the RluA family, showed a slight upregulation upon PUS1 knockout (Supplementary Fig. 3g). We noted that in PUS1-knockout cells, PUS1 mRNA was not lost despite the loss of protein function, likely because the premature stop mutation did not fully trigger nonsense-mediated RNA decay. Following PUS1 knockout, the most prominent increases in Ψ modification were observed at Ψ30 and Ψ31. tRNAs exhibiting PUS1-dependent increases in Ψ at U30/31 contained a distinct “RGU (R=G or A)” motif, whereas non-interplay tRNAs harbored a “CKU (K=G or U)” motif, indicating that distinct sequence features might account for their differential susceptibility to PUS1-mediated interplay (Supplementary Fig. 3h). Although Ψ31 and Ψ32 have been reported to be catalyzed by RluA family enzymes in yeast40,41, the substrate specificity of the RluA family in humans remains unclear. To further investigate the molecular basis of this interplay, we subsequently analyzed the potential tRNA substrates of human RluA family members.

Characterization of tRNA targets for the human RluA family enzymes
The human RluA family consists of four proteins: RPUSD1, RPUSD2, RPUSD3, and RPUSD4. It has been reported that RPUSD4 mediates pseudouridylation at Ψ3067 of mitochondrial 16S rRNA42; however, the substrates of the other three proteins remain uncharacterized. We generated RPUSD1−/−, RPUSD2−/−, and RPUSD3−/− knockout cell lines using CRISPR/Cas9 in HEK293T cells. The knockout efficiency was confirmed by genomic DNA sequencing, demonstrating efficiencies of ~100.0% for all three cell lines (Supplementary Fig. 4).

We then utilized PRAISE to identify the modification sites and levels of Ψ targets by comparing Ψ sites of wild-type and RPUSD KO cells. Using the stringent criteria (modification level difference >5.0% and reduction in the modification level of Ψ sites >2-fold), we identified 16 RPUSD1-dependent Ψ sites in tRNA-Arg and tRNA-Ile (Fig. 4a and Supplementary Data 7). In contrast, we did not find any confident cytosolic tRNA substrates for RPUSD2 and RPUSD3 proteins (Fig. 4b, c). In mitochondrial tRNA, we identified 3 RPUSD2-dependent Ψ sites, and the modification levels of these sites were nearly abolished (dropped to 0.00%) after RPUSD2 knockout (Fig. 4e and Supplementary Data 7). In contrast, neither RPUSD1 nor RPUSD3 was found to modify Ψ targets in mitochondrial tRNA (Fig. 4d, f). Given that the subcellular localization of PUS proteins can influence their accessibility to potential RNA targets, we used MitoFates43 to predict the mitochondrial localization of the RPUSD proteins. Our results demonstrated mitochondrial targeting signals were absent in the N-terminal of RPUSD1 isoforms (Supplementary Fig. 5a), whereas they were indeed present in the N-terminal of two RPUSD2 isoforms (Supplementary Fig. 5b), suggesting that RPUSD2 is localized in mitochondria. Collectively, RPUSD1 was found to modify the cytosolic Ψ targets and RPUSD2 was found to modify the mitochondrial Ψ targets.

We further analyzed the modification sites and levels of these Ψ targets. Interestingly, in addition to modifying Ψ30 of cyto-tRNAs, RPUSD1 also modified Ψ72 of cyto-tRNA (Fig. 4g). Specifically, the modification level of Ψ30 on cytosolic tRNA-Ile-AAT decreased from 20.0 to 0.0% after RPUSD1 knockout, and the modification level of Ψ72 on cytosolic Arg-TCG tRNA dropped from 58.1 to 0.0% after RPUSD1 knockout (Fig. 4h and Supplementary Fig. 6a). Consequently, we unexpectedly identified Ψ72 as a site specifically catalyzed by RPUSD1, which is entirely distinct from any RluA family homolog in yeast40,41. RPUSD2-dependent Ψ sites are located at U31 and U32 of mitochondrial tRNA (Fig. 4i). The enzyme responsible for Ψ31 modification in mitochondrial tRNAs is highly similar to its yeast homolog Pus940. Notably, our results revealed that RPUSD3 has no targets on any cyto-tRNAs or mt-tRNAs, indicating it lacks catalytic activity due to the absence of a critical aspartate residue2. Furthermore, we investigated whether these knockouts or targeted Ψ sites affected tRNA expression. Analysis revealed that there is no significant impact on the global expression levels of either cyto- or mt-tRNAs (Supplementary Fig. 6b, c). The impact of Ψ targets modified by the human RluA family in tRNA biogenesis and translational fidelity warrants further investigation.

Interplay of RPUSD1 and PUS1
The results showed that RPUSD3 does not act on substrates in either cytosolic or mitochondrial tRNAs; thus, the observed increase in Ψ30 and Ψ31 sites following PUS1 knockout was not caused by altered RPUSD3 expression. Interestingly, upon RPUSD1 knockout, in addition to a significant reduction in Ψ30 on its tRNA substrate Ile-AAT, we observed a notable increase in the adjacent modification site, Ψ28, on the same tRNA molecule (Fig. 4j). Specifically, our quantitative data showed that the Ψ level of Ψ28 increased by more than 30% following RPUSD1 knockout. To determine whether this effect was due to upregulated PUS1 expression leading to increased modification levels of its targets, we examined PUS1 expression levels. However, PUS1 expression remained unchanged after RPUSD1 knockout (Supplementary Fig. 6d), indicating that the observed changes in modification levels were not attributable to altered PUS1 expression. Given the interplay between adjacent Ψ sites, we speculate that this phenomenon may arise from steric hindrance among different PUS enzymes and may reflect a broader regulatory mechanism involving multiple Ψ-modifying enzymes within the anticodon stem-loop.

To assess whether PUS interplay depends on the catalytic activity of PUS1, we performed rescue experiments in PUS1 knockout cells using wild-type PUS1 or the catalytic-dead D118A mutant11 and the disease-associated nonsense variant E192X10 (Supplementary Fig. 7a–e). Specifically, PUS1 knockout HEK293T cells were transfected with constructs encoding wild-type human PUS1 isoform 2 or the corresponding D118A or E192X mutant proteins. PRAISE sequencing showed that wild-type PUS1 successfully restored the Ψ28 modification and reduced the elevated interplay-associated Ψ30 level in Ile-AAT tRNA. In contrast, neither D118A nor E192X rescued Ψ28 modification or the interaction with Ψ30 (Supplementary Fig. 7f–h). These results indicate that loss of catalytic activity, whether caused by mutation of the catalytic residue or by pathogenic truncation, might abolish the ability of PUS1 to mediate PUS interplay.

The interplay of PUS7 with PUS1 and PUS3
Since Ψ35 in tRNA has been reported as a PUS7 target44, we then aimed to explore potential interplays between tRNA modification enzymes PUS1 and PUS7. To this end, we applied PRAISE to the PUS7−/− cell line to identify its substrate sites. We first found that PUS7 knockout does not affect the expression of other PUS enzymes (Supplementary Fig. 8a), nor does it alter the global expression levels of cytosolic or mitochondrial tRNAs (Supplementary Fig. 8b, c). Then we identified 100 PUS7-dependent Ψ sites in cytosolic tRNAs (Fig. 5a and Supplementary Data 1). Their modification levels decreased significantly by >50%, and their contexts were consistent with the well-known PUS7 motif “UGUAG” (Fig. 5b, c). We also found that PUS7 mainly mediated the formation of Ψ13, Ψ20b, and Ψ50 in 83 cytosolic tRNAs (Fig. 5d). Additionally, we did not identify any Ψ modification sites in mt-tRNAs, indicating the absence of PUS7’s localization and catalytic activity in the mitochondria.

After PUS7 knockout, we observed that the levels of Ψ28 in Arg-TCT tRNAs significantly increased from 55.4 to 70.1%, while its target site Ψ36 decreased significantly from 11.2 to 0.0% (Fig. 5e). Conversely, in PUS1−/− cell line, the modification levels at Ψ28 sites decreased from 55.4 to 10.8%, while the modification levels at Ψ36 sites increased from 11.2 to 56.6%. These results provide evidence for the interplay between PUS1 and PUS7 in the anticodon stem-loop region. We then analyzed the changes in the levels of nearby targeted Ψ sites and found that the modification levels of Ψ35 sites in Tyr-GTA tRNA isoacceptors decreased from 87.7 to 0.0%, while the modification levels of Ψ39 sites increased from 13.0 to 62.8% (Fig. 5d, f). Since Ψ39 is located in the anticodon stem-loop and has been shown to be catalyzed by PUS345, these results suggest an interplay between PUS7 and PUS3 in the anticodon loop of cyto-tRNAs. tRNAs carrying U39 were grouped into PUS1- or PUS7-targeted and non-targeted sets. In both PUS1 and PUS7 knockout cells, Ψ39 increased only in the respective targeted tRNAs, while non-targeted tRNAs showed no significant change (Supplementary Fig. 8d–g). These results indicated that the observed Ψ39 elevation occurs specifically at PUS-dependent substrate sites rather than through a global increase effect.

To further investigate whether the interplay among PUSs occurs at other positions of cyto-tRNA, we applied PRAISE to the generated TRUB1−/− cell line, where TRUB1 catalyzes the highly-modified Ψ55 sites in the TΨC loop of tRNA. We first detected that TRUB1 knockout did not affect the expression of other PUSs (Supplementary Fig. 9a). Using the criteria mentioned above, a total of 28 TRUB1-dependent sites were identified (Supplementary Fig. 9b and Supplementary Data 1). These sites were primarily localized at Ψ55, and their context conformed to the TRUB1 motif “GUUCNA” (Supplementary Fig. 9c, d). Compared to PUS1-dependent and PUS7-dependent sites, the reduction in modification levels of TRUB1-dependent sites was less pronounced. To investigate the potential redundancy between TRUB1 and PUS10 in cytosolic tRNAs46,47, we utilized a previously generated PUS10−/− HEK293T cell line21 and performed PRAISE analysis on small RNAs. Consistent with our earlier CMC-based findings, Ψ54 sites in tRNA were completely abolished upon PUS10 knockout. Moreover, 21 Ψ55 sites showed substantial decreases in the PUS10-/- cells, including those in Ala-AGC, Val-AAC, Val-TAC, and Val-CAC tRNAs (Supplementary Fig. 9e). Notably, we also identified an additional set of 14 Ψ55 sites that exhibited modest reductions and were affected by both PUS10 and TRUB1 knockout (Supplementary Fig. 9f). Together, these results revealed functional redundancy between PUS10 and TRUB1 at cytosolic tRNA Ψ55 sites.

In contrast to cytosolic tRNA, we identified 6 TRUB1-dependent sites, all located at Ψ55, in mitochondrial tRNAs (Supplementary Fig. 9g and Supplementary Data 1), and their modification levels were eliminated upon TRUB1 knockout (Supplementary Fig. 9h). These findings excluded Ψ55 synthase activity of TRUB247. Consistent with its mitochondrial role48, immunofluorescence staining in HeLa cells showed that TRUB2 is exclusively localized to mitochondria (Supplementary Fig. 10a), raising the possibility that TRUB2 may contribute to mitochondrial function via non-catalytic mechanisms. We found that TRUB1 knockout did not affect the expression of cyto- and mt-tRNAs, which confirmed that the change of Ψ level was not caused by tRNA expression (Supplementary Fig. 9i, j). Meanwhile, although the modification level of Ψ55 on cyto-tRNA and mt-tRNA decreased, the modification level of nearby sites did not significantly change (Supplementary Fig. 9d). Therefore, the results suggested that the interplay among different PUSs in tRNA primarily occurs in the anticodon stem-loop region.

tRNA substrates for PUS1L and PUS7L
To further expand the tRNA Ψ landscape, we generated knockout HEK293T cells for stand-alone PUS1L and PUS7L and performed PRAISE sequencing to identify their Ψ substrates. Remarkably, PUS1L knockout eliminated all seven Ψ39/40 sites in mitochondrial tRNAs while not significantly influencing Ψ sites in cyto-tRNAs, indicating that PUS1L exclusively modifies mt-tRNAs (Supplementary Fig. 10b and Supplementary Data 1 and 2). Immunofluorescence analysis also showed its mitochondrial localization (Supplementary Fig. 10c, d). In contrast, PUS7L knockout led to a pronounced reduction in 25 Ψ sites (Supplementary Fig. 10e, f). These PUS7L-dependent sites were restricted to type II tRNAs, which are defined by an extended variable loop. The affected tRNAs included Ser-GCT, Leu-AAG, Leu-CAA, and Leu-TAG, and all Ψ sites were located within the variable loop region (Supplementary Data 1). Together, these results indicated that PUS1L and PUS7L display substrate specificity entirely distinct from their paralogs, PUS1 and PUS7, suggesting that they may have specialized functions that warrant further investigation.

Methods
Cell culture
HEK293T cells were maintained in DMEM medium (Corning) supplemented with 10% FBS (Gibco), 1% GlutaMAX (Gibco), and 0.5% penicillin-streptomycin (Gibco) at 37°C with 5% CO2. Cells within passage 6 were used for experiments. Stress conditions were applied as follows: heat shock was induced by incubating cells at 43°C for 2 h, followed by a 6 h recovery at 37°C, and oxidative stress was induced by treating cells with medium containing 1 mM H₂O₂ (Sigma) for 1 h.

Generation of stable knockdown and CRISPR knockout cell lines
The shRNA targeting DKC1 (TRCN0000010325) was designed and transfected into HEK293T cells as described in our previous study23. PUS1, PUS7, TRUB1, PUS10 KO in HEK293T cells were generated by previous studies17,21,23,61. RPUSD1, RPUSD2, RPUSD3, PUS1L, and PUS7L KO HEK293T cells were newly generated by CRISPR–Cas9 technology62. sgRNA sequences were listed as follows:
RPUSD1-sgRNA: CCGGCAGCGCGTACAGGTGC
RPUSD2-sgRNA: AAGGTCGTAGCGCCAATGTC
RPUSD3-sgRNA: CTCCAGAACCGGCCCAAAAC
PUS1L-sgRNA: AGGCCCTCAACACACACCTG
PUS7L-sgRNA: CTGAACCAGACTGATGATTT

Plasmid construction and cell transfection in rescue experiments
Wild-type human PUS1 isoform 2 was cloned from cDNA into the pcDNA3.1(+) vector (Invitrogen). Using this plasmid as a template, we generated two mutant variants: a premature termination mutant10 (PUS1-iso2 E192X) and a catalytically inactive mutant11 (PUS1-iso2 D118A). Site-directed mutagenesis was performed using the Mut Express II Fast Mutagenesis Kit V2 (Vazyme) with the following primers. Plasmids were transfected into PUS1-KO HEK293T cells using Lipofectamine LTX with Plus™ Reagent (Invitrogen). Cells were harvested 48 h post-transfection for RNA and protein extraction.
PUS1-iso2 E192X_F primer: TCAGGTTAGACCTACCGCCTGAGCGCCGAGA
PUS1-iso2 E192X_R primer: GGTAGGTCTAATCCTGAACGTCCCGGTCCTTG
PUS1-iso2 D118A_F primer: GACAGCCAAGGGTGTGTCCGCAGCCGGCCAGG
PUS1-iso2 D118A_R primer: ACACACCCTTGGCTGTCCGGGCGCAGCGCTGG

Small RNA purification
Total RNA was extracted with TRIzol (Life technologies) followed by isopropanol precipitation, according to the manufacturer’s instructions (Invitrogen). The resulting total RNA was treated with DNase I (NEB) to avoid DNA contamination. Size selection was performed using MEGAclear Transcription Clean-Up Kit (Ambion) to deplete large RNA. The small RNA fraction was subsequently purified by ethanol precipitation.

PRAISE for small RNA
About 300 ng of small RNA from each technical replicate was denatured and demethylated by incubating with 2× molar ratio of wild-type AlkB, and subsequently recovered through phenol-chloroform extraction and ethanol precipitation. About 100 ng demethylated small RNA was transferred into a new nuclease-free tube and used as “untreated” and stored at −80°C. Sulfite/bisulfite reagent was freshly prepared by dissolving 2.58 g of potassium sulfite (Sigma-Aldrich) and 0.30 g of sodium bisulfite (Sigma-Aldrich) in 8 ml RNase-free water. The 100 mM hydroquinone was prepared freshly by adding 11.01 mg of hydroquinone (Sigma-Aldrich) to 1 ml of RNase-free water. For sulfite/bisulfite treatment, 5 µl of RNA was mixed with 50 µl of sulfite reaction, which is a 100:1 mixture of sulfite/bisulfite reagent and 100 mM hydroquinone, and subjected to heat incubation at 70°C for 5 h. The reaction mixture was desalted by twice passing through Micro Bio-spin 6 chromatography columns (Bio-Rad). Desalted RNA was then incubated with an equal volume of 1 M Tris-HCl (pH 9.0) at 75°C for 30 min. The reaction was then immediately stopped by chilling on ice and ethanol precipitation.

After ethanol precipitation, “treated” and “untreated” small RNA samples were subjected to library construction according to the eCLIP protocol with several modifications63. RNA samples were dephosphorylated with T4 PNK (NEB). 3' adapter (5'-rApp-AGATCGGAAGAGCGTCGTG-C3spacer-3') was ligated using T4 RNA Ligase 2, truncated KQ (NEB) at 25°C for 2 h. The excess adapter was digested using 5' deadenylase and RecJf (NEB). RNA was reverse transcribed using Maxima H minus (Invitrogen) and RT primer (ACACGACGCTCTTCCGATCT), following by RNase H treatment. 5' linker (5Phos-NNNNNNNNNNAGATCGGAAGAGCACACGTCTG-3SpC3) ligation was performed with T4 RNA Ligase 1, High Concentration (NEB). The ligated cDNA was purified using silane beads and then subjected to PCR amplification using PCR primers: 5'-AATGATACGGCGACCACCGAGATCTACACTCTTTCCCTACACGACGCTCTTCCGATCT-3', 5'-CAAGCAGAAGACGGCATACGAGATXXXXXXGTGACTGGAGTTCAGACGTGTGCTCTTCCGATC-3' (XXXXXX represents index sequence). PCR products were purified using 8% native-PAGE gel. The libraries were sequenced on Illumina HiSeq X10 with a paired-end 2×150 bp read length.

Preprocessing of raw reads
The PRAISE library construction enables reads 2 to retain the original RNA direction, so only reads 2 was used for analysis in this study. First, Cutadapt (version 4.8)64 was used to remove the adapter and perform quality control of the sequencing reads. The key parameters are: cutadapt -j 28 --times 1 -e 0.1 -O 3 --quality-cutoff 25 -m 55. After that, based on the 10 bp unique molecular identifier (UMI) at the 5' end of reads 2, the duplicated sequences were removed by Seqkit (version 0.13.2)65. The key parameters are: seqkit rmdup -s. Finally, the 10 bp UMI was extracted by umi_tools (version 1.0.0)66. The key parameters are: umi_tools extract --extract-method=string --bc-pattern=NNNNNNNNNN. In addition, the library construction added 6 bp bases at the 3' end of insertion, which was also removed by umi_tools (version 1.0.0). The key parameters are: umi_tools extract --extract-method=string --3prime --bc-pattern=NNNNNN.

Cleaned reads alignment of PRAISE
Human cytosolic tRNAs reference (GRCh38) was retrieved from the database GtRNAdb67. Only non-redundant tRNA sequences were retained, and a “CCA” sequence was added to the 3' end of each tRNA. The mitochondrial tRNAs reference was downloaded from NCBI. Other small RNA references (snRNA, snoRNA, microRNA, yRNA, and <200 nt rRNAs) originated from the human transcriptome (GRCh38). The snRNA and snoRNA entries were streamlined to only one representative sequence per type.

A two-step strategy was implemented to ensure the alignment accuracy: (1) BWA-MEM was used to align the cleaned reads to the reference transcriptome of total small RNAs. The key parameter is: bwa mem -Y. (2) A custom Python program was used to realign the mapped reads from the previous step. This program was based on biopython and adopted an optimized scoring matrix to precisely resolve deletion signals. The key parameters were: python realignment.py –fast -ms 4.8. After that, the custom Python program remove_end_signal.py was used to remove 3'-end deletion signals to reduce false-positive sites. Samtools (v1.14) was used to convert the BAM file into an mpileup format file68. The key parameters are: samtools mpileup -d 1500000 -BQ 0 --ff UNMAP,QCFAIL -aa. Finally, the information in the mpileup file was counted and converted into a BMAT format file for subsequent signal identification.

Identification of Ψ signal
Ψ site identification required three criteria: (1) For each T site, a 2×2 contingency table was constructed using deletion counts and non-deletion counts (total counts - deletion counts) from treated and untreated samples. Significant inter-group differences were assessed by χ² or Fisher’s exact test. The p value is required to be <0.0001; (2) The deletion rate was defined as follows: deletion counts/total counts. The deletion rate difference between treated and untreated samples is required to be >5%. (3) Ψ sites must be detected in both technical replicates. The mean detected signal between replicates was reported as the final signal of each Ψ site.

Identification of PUS-dependent Ψ sites in small RNA
To identify PUS-dependent sites, the final Ψ list in WT samples was taken as a reference. The decreasing rate of deletion signal in tRNA between PUS-KO samples and WT samples is required to be >50%. The decreasing rates of deletion signal in snRNA and snoRNA between DKC1-KD samples and WT samples are required to be >30%. The total coverage of PUS-dependent sites is required to be >10 reads.

Analysis of tRNA expression
Adapter trimming and quality control of raw reads were performed using Cutadapt (version 4.8). The cleaned reads were mapped to the tRNA reference (GRCh38) using BWA-MEM with default parameters. The mapped reads were processed using samtools (v1.14). RNA coverage was quantified using pysam. RNAs with low technical reproducibility were discarded. The criteria is: log2 fold change (coverage of technical replicate 1/coverage of technical replicate 2) >1. tRNA expression levels were calculated using RPM.

Analysis of RNA expression
Cutadapt (version 4.8) was used for quality control and adapter trimming of raw reads. The cleaned reads were mapped to the reference transcriptome (GRChg38) using Hisat2 with default parameters69. Mapped reads were processed by samtools (v1.14). Gene coverage was quantified by featureCounts. The key parameters are: featureCounts -t exon -g gene_id. RNAs with low technical reproducibility were discarded. The criteria is: log2 fold change (coverage of technical replicate 1/coverage of technical replicate 2) >1. The expression of PUSs was calculated by FPKM.

Motif analysis, correlation analysis, and statistics analysis
We extracted the upstream and downstream sequences of the target site, and then ggseqlogo in R (v. 4.3.1) was used for Motif visualization. The data in this study were continuous data, so the statistical parameter in the correlation analysis was the Pearson correlation coefficient. The statistical tests used in this study mainly included the two-sided Student’s t-test and the Wilcox rank-sum test. The two-sided Student’s t-test was used when the data conformed to the normal distribution and the data variance was homogeneous; otherwise, the Wilcox nonparametric rank-sum test was used. Error bars represent mean ± s.d.

Target prediction on snRNA and snoRNA
We performed an in-silico target prediction analysis based on the snoRNA database36. For each pseudouridine site detected, a 21-nt window containing 10 nt upstream and downstream of the modified position was extracted to generate the target set. The 3' and 5' catalytic hairpin regions of all H/ACA box snoRNAs were obtained using RNAfold (2.7.0) from the ViennaRNA Package70 and used as the snoRNA query library. We then used RNAsnoop (2.7.0) from the ViennaRNA Package to generate specific binding models of the query library and target set. The key parameters are: RNAsnoop -s snoRNA query library -t target set -d 2. We employed filtering criteria in RNAsnoop to minimize false positives: (1) ≥4 base pairs in the 5' duplex, (2) ≥4 base pairs in the 3' duplex, and (3) a total predicted free energy <–10 kcal/mol.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "This study comprehensively profiled pseudouridine (Ψ) modifications in small RNAs—including cytosolic and mitochondrial\ntRNAs, snRNAs, and snoRNAs—in HEK293T cells using the quantitative Ψ detection method known as PRAISE. As a\nresult, the authors newly identified the RNA substrates and corresponding modification sites for multiple pseudouridine\nsynthases (PUSs). Furthermore, they demonstrated crosstalk among multiple PUS enzymes in the anticodon stem-loop\nregion of tRNAs, whereby these enzymes influence each other's modification activity. These findings represent the first\ncomprehensive, quantitative, and high-resolution dataset of Ψ modifications and mark a significant milestone in the field of\nRNA modification research. Below are several points for consideration.\nAdding a discussion on the physiological relevance of the observed interplay among PUS enzymes would greatly enhance\nthe manuscript. For instance, could such interplay reflect a compensatory mechanism that maintains tRNA functionality in\nresponse to the loss of specific modifications?\nIn PUS1 and PUS7 knockout cells, increased Ψ modification at position 39 was observed in many tRNAs that contain the\nrespective target sites for these enzymes. Do similar changes also occur in tRNAs that do not carry these target sites?\nDiscussing this point would be important for understanding potential interactions among PUS enzymes.\nTo provide a more comprehensive understanding, it would be valuable to present data on the Ψ modification sites and their\nlevels in snRNAs and snoRNAs in RPUSD1–3 knockout cells.\nIt is useful to summarize all the data (Ψ position, tRNA species, and dependency on each enzyme) in a single table (e.g., an\nExcel sheet).\nOther comments\nThe data indicate that Ψ8 is present only in Gly-GCC, but Reference 14 claims that Ψ8 is found in fragments derived from\nAla, Cys, and Val as well. How do the authors reconcile this discrepancy?\nFor abbreviations such as CMC, the full name must be given on first appearance.\nLine 118: Does this refer only to cytosolic tRNAs? Please clarify whether mitochondrial tRNAs are included.\nLine 159: Describe the specific values (total number of sites, etc.) on each region.\nLine 162-163~. The number of Ψ sites at each position is described, but the meaning of this number is difficult. Because Ψ is\ninduced only in uridine, this number is significantly affected by the ratio of uridine at each position. Describing the uridine\nnumber at each position is helpful.\nLine 171: Please clarify this ‘median level’ is calculated from all tRNAs or modified tRNAs only.\nLine 242: Instead of “UGG”, would “USG (S = G or C)” or “USR (R = G or A)” be more appropriate?\nLine 253: “Fig. 2c,d” should be corrected to “Fig. 2b,c”.\nLine 481: The sentence appears to be unclear or grammatically incorrect. Please revise for clarity.\nLine 482: This line is missing or cut off, making it unreadable.\nLine 488: Please provide the DNA sequences used for PCR.\nLines 497 and 499: Should “8 bp” be corrected to “10 bp”?\nLine 540: Please clarify what “10” refers to. 10 reads?\nFigure 1a: Are the values the average of two biological replicates?\nFigure 1b: What are the different sites between different repeats? Show the details.\nFigure 1e: The label “arm and loop” should be corrected to “stem loop” as described in the text.\nFigure 1e legend: Please provide an explanation of the colored bars.\nFigure 1e, 3e, 5d: The heat maps are very hard to read, not only X axis, Y axis and signal, please revise them to a more\nclear version.\nFigure 1f: “ 8” is not shown. It would also be helpful to indicate secondary structure information (e.g., anticodon stem loop) to\ncorrespond with Fig. 1e.\nFigure 2: Please define the abbreviations “UT” and “NC” in the figure legend. In addition, ensure consistent labeling—for\nexample, “DKC1-NC cells” is used in the legend, whereas “DKC1-NC-UT” appears in the figure.\nFigure 3e,5d: Overall color of this figure is a light pink. Is that OK?\nFigure 5g: Is it correct that positions 27 and 35 of Arg-TCT should actually be 28 and 36? Similarly, should position 28 of\nMet-CAT be 27?\nEx. Fig. 1a: “Alkb” should be corrected to “AlkB”.\nEx. Fig. 1b: Should “puro” be replaced with “spacer”?\nEx. Fig. 1b: Should “PNK” be replaced with “CIAP”?\nEx. Fig. 1e: It would be helpful to include a correlation plot between PRAISE and BACS.\nEx. Fig. 1f: A figure legend is missing.\nEx. Fig. 1f–g: Please correct the label to “1g–h”.\nEx. Fig. 2a and 2c: Clearly indicate whether the data show cytosolic or mitochondrial tRNAs. The figure legend should\nspecify “left” or “right” accordingly. The same applies to Ex. Fig. 5b and 5c, and Ex. Fig. 6b and 6c.\nEx. Fig. 2d, Extended Fig. 5d, Extended Fig. 6a, and Extended Fig. 7a: In these figures, it is not immediately clear that the\ntargeted gene itself is not lost in the knockout condition. It would be helpful to add an explanatory note in the main text.\nEx. Fig. 3: The DNA sequence is shown, but it is unclear which region of the gene it represents. Please indicate where the\nindel occurred in the gene structure. Also, figure legend labels (a–f) do not match the figure panels. Was Western blotting\nperformed?", "reviewer_2": "In this manuscript entitled “Quantitative analysis of small RNA pseudouridylation reveals interplay of PUS enzymes in tRNA\nanticodon stem-loop”, Liu et al adapted their original PRAISE method to small RNA Ψ profiling by incorporating AlkB\ndemethylation to remove the methylation of several tRNA modifications to improve the readthrough in RT reaction, enabling\nquantitative and site-specific profiling of Ψs in small RNA. The method worked very well for quantitative profiling the Ψ sites\nin small RNA including tRNA, snRNA and snoRNA. By perturbing seven PUS enzymes using stable knockdown or\nknockout cell lines, the authors made interesting discoveries including: (1) pseudouridylation of snoRNAs can be mediated\nother stand-alone enzyme PUS7 in addition to DKC1; (2) quantification of the Ψ fractions in KO and WT samples not only\nallow the assignment of the tRNA substrates of individual PUS enzyme, but also reveal unexpected interplay among PUS\nenzymes in shaping Ψ patterns within the ASL; and (3) for the first time the tRNA targets of the human RluA family were\ndetermined the human tRNA targets which differs from that of its yeast counterparts. The manuscript was written very clearly\nand the evidence they showed are adequate to support their conclusions. Their finding will definitely benefit more functional\nstudies of Ψ and PUS enzymes. Therefore, this reviewer supports its publication in Nature Communications after minor\nrevision.\nThe following are some minor comments:\nLine 107-108, “However, specific modifications, including m1A and m3C\nmodifications, can block reverse transcription”. Besides m1A and m3C, how about other modifications in tRNA such as m1G\nand m22G? Can these modifications also block RT reactions? Can these RNA modifications be demethylated by AlkB\ntreatment?\nLine 118-9: “rRNA accounted for 12.3%”. The small RNA fraction should not have much rRNA contamination, is there some\nrRNA degradation in the course of small RNA purifications?\nLine 123 “In the cyto-tRNAs, 713 and 714 Ψ sites were detected in two technical replicates, respectively, with 705\noverlapping Ψ sites, thus demonstrating high reproducibility…” it would be nice to describe the quality of these libraries,\nsuch as the mapping ratio, deduplication rate, insert sizes, and the criteria to call these sites. In addition, for the non-\noverlapped Ψ sites, are their fractions lower in general than the detected sites? What might be the reasons why these sites\nfailed to overlap?\nLine 174-175 “modification level, reaching a median of 83.4%. Ψ35, located at the center of tRNA-Tyr anticodon, showed a\nhigh modification level of 86.9%”. The detected modification level is calculated based on the deletion rate and the calibration\ncurves?\nLine 213 “12 high-confidence DKC1-dependent sites among 28 Ψ sites in snRNA”, for the other 16 Ψ sites in snRNA, any of\nthem is sensitive to the KO of the other PUS enzymes?\nLine 411-414, “facilitate characterization of additional PUS enzymes, such as PUS1L and PUS7L, enabling comprehensive\nmapping of tRNA Ψ modifications and advancing our understanding of PUS-mediated Ψ formation.” Is there any specific\nreason why the authors don’t include the KO of PUS1L and PUS7L in this study?", "reviewer_3": "Major comments:\n- The authors performed a comprehensive analysis of pseudouridylation in human tRNAs as well as other small RNAs\nincluding snoRNAs and snRNAs. The data demonstrate all of the positions of pseudouridylation, while quantifying levels\nusing their previously established method, named PRAISE. They were also able to assign many of the sites to specific\npseudouridine synthases through analysis of cells with knockdown or knockouts of 7 different pseudouridine synthases. The\ninnovation derives primarily from the large scale of data, and this work will serve as a very useful resource for the community\nstudying RNA modifications in tRNAs, snRNAs, and snoRNAs. The authors do highlight some novel findings, including a\nsnoRNA pseudouridine site not catalyzed by DKC1 but by “stand-alone” synthase PUS7, as well as new examples of\napparent interplay between different pseudouridines in the anticodon stem/loop, in which loss of one site can lead to either\nincreases or decreases of other sites catalyzed by other enzymes. These examples are limited in comparison to the rest of\nthe data presented, however, they are not orthogonally validated by another method and are not considered in regard to\nhypotheses or models that could be tested, which would have extended their significance for the field. This work represents\na great resource for the field, but alone does not present a significant new understanding of tRNA biology, although the data\nhave this potential if extended through further analysis or discussion.\n- This work showcases the potential of PRAISE to establish the record of where pseudouridine is and at what levels\nquantitatively. A limitation is that information about other modifications is not recorded. This may not influence interpretation\nof analysis of snoRNAs and snRNAs, because those RNAs are not known to contain a wide variety of modifications, at high\ndensity. But for tRNAs, there is a lot of lost information in not considering how loss of pseudouridine may influence other\nmodifications. For example, several examples are known in which loss of pseudouridine in tRNA at specific sites can affect\nother modifications (Porat 2025 for review). While it remains challenging to record this information, it is suggested that the\nauthors nonetheless acknowledge this limitation of the approach in the text, and discuss the possibility that some of the\ninterplay that they observe is “indirect” (in that it is happening through one pseudouridine influencing another modification,\nwhich is in turn influencing a different pseudouridine site). To summarize, due to the complexity of the tRNA modification\nlandscape, it is difficult to draw conclusions for interplay between modifications when only examining one type out of more\nthan two dozen possibilities in tRNA.\n- The quantitative analysis of pseudouridylation site modification is impressive, but sometimes other aspects of quantification\nneed improvement. For example, in the analysis presented in lines 227–236, it describes how 48/73 sites are DKC1-\ndependent. Then after analysis of three mutants, it is concluded that one site is PUS7-dependent. Was this site one of the 25\nthat were not DKC1-dependent? Presumably? But what about the other 24? Are there any ideas for what those could be\ncatalyzed by? Given that this result was stated in the abstract, there needs to be a more quantitative description of the\nnumber of sites in snoRNAs that are DKC1-dependent, and the number of new sites that have been revealed to be\ndepended on a stand-alone synthetase, which from this study appears to be one (SNORD65). Moreover, because of the\napparent significance of this result—it is mentioned in the abstract—the authors can consider whether they would like to\ndemonstrate it by an orthogonal approach. Given that the entire study basically consists of PRAISE assays, this type of\nresult is of the type that may most benefit from an orthogonal validation.\n- The description of changes on mt-tRNAs resulting from the absence of a TRUB1 were interesting but would benefit from\nthe addition of a description of the possible role of TRUB2. It would also be beneficial to know whether the cytosolic\ntRNAs that were not changed by TRUB1 deletion lacked the motif described in line 275. Or differed in other features. The fact that a\nrelatively small number of these isoacceptors were affected, and to incomplete levels is interesting, and more discussion of\nPUS10 vs. TRUB1 would improve impact of understanding these interesting these results.\n- In the case where there is “interplay” between different modified sites in the anticodon loop, how is “labeling efficiency” of a\ngiven site affected when there is a strong decrease in a neighboring site? Presentation of data that would refute such an\nartifactual model for changes would be appreciated.\n- Pathogenic variants of PUS1 and PUS7 are included in the introduction, but not in the discussion. Given that this study\nfound several examples where PUS1 or PUS7 catalyzed pseudouridine sites impact other sites, it would be reasonable to\nspeculate how these interactions may be involved in a disease mechanism. In the case of PUS1 disease variants, several\nare known to have a reduction/loss of catalytic activity, which could possibly cause altered interactions with other\npseudouridine sites compared to a PUS1 knockout. In summary, connecting the findings of this study to pathogenic variants\nwould emphasize the broader importance of the results.\n- The paper title and results section are centered on tRNA. The section on snoRNA and snRNA, while adding value to the\nstudy, feels out of place. It is suggested to move this section to the end to maintain the focus on tRNA that appears to be the\nintent of the authors. And to reframe the study as a focus on tRNA, but that additional analyses of these other classes of RNA\nadded further novelty to the study.\n- The analysis of mRNA and tRNA expression are not sufficiently explained. Why were these done? How do they facilitate\ninterpretation of the PRAISE data? Some examples from the literature should be provided that defend the idea that loss of\none PUS enzyme or other RNA modification enzyme has been shown to lead to changes in another enzyme, and to a\nsignificant enough degree that it would affect the activity of that second enzyme. Otherwise, this analysis seems to yield\nlimited value for the space that it occupies, and without further discussion to emphasize its importance.\n- A major focus of this study is the interplay between pseudouridine sites in the anticodon stem-loop, of which they found\nseveral interesting examples. However, there are many tRNAs that are modified by PUS1 and do not have an increase in Ψ\nlevels at U30/31 (Fig. 3e). Are there any shared features (i.e. tRNA sequence, structure, etc.) that may influence whether an\ninteraction between Ψ27/28 and Ψ30/31 is present? In lines 386-387, the authors state that interactions between PUSs\nprimarily occurs in the anticodon stem-loop, which is somewhat contradictory to the data presented in Fig. 3e and 5d that\nshow several examples of interactions between Ψ sites in other regions of the tRNA. The mechanisms proposed in lines\n429-439 could benefit from additional references to structural/biochemical literature on PUSs.\n- Some of the figure legends are very short and do not provide the complete information necessary to interpret the plot. For\nexample, there is insufficient explanation of the heat maps shown in 1e, 3e and 5d/f. For example, it is not clear the order\nthat the isodecoders are labeled, and in general there just needs to be more detail in how these are presented.\nMinor comments:\n- Figure 3F, 4J: the labels caused some confusion. Is this a standard nomenclature for listing a tRNA isodecoder and its\ncorresponding position? Simply adding the word “position” prior to the last number would distinguish the gene name and the\nposition in the tRNA of those mapped reads. Alternatively, the data can be split by position, for example using brackets on\nthe right side of the plot to group all “position 27” rows together.\n- line 427–429: these citations are incomplete in presenting what is known about how psi55 influences m1A58. There are\nnumerous cases known in yeast in which the influence is positive and not negative (Barraud et al 2019, Shaw and Thomas\net al 2024, reviewed in Porat 2025).\n- lines 251–252: the meaning of this statement is not clear. Why would a difference in tRNA expression between two\nsamples cause difficulties in assignment of the modified positions in tRNA?\n- Figure 5D: The figure legend describes the plot as a histogram, when 5D is a heatmap. Descriptions for 5D and 5E may\nhave been swapped in legend.\n- Figure 5G: this illustration and legend need to be revised for clarification. There are no arrows being shown, and the\nidentity of starting and terminal nodes is not clear and thus direction of change. The first two simply have solid lines\nconnecting the two sites.\nThe Github link does not work and the GEO accession number states that the data will not be public until 2027. Therefore, it\nis not possible to reproduce the analyses with these current circumstances.", "reviewer_4": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature\nCommunications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career\nResearchers who co-review manuscripts."}, "author_response": "Point-by-point Response to Reviewer #1\nThis study comprehensively profiled pseudouridine (Ψ) modifications in small\nRNAs—including cytosolic and mitochondrial tRNAs, snRNAs, and snoRNAs—\nin HEK293T cells using the quantitative Ψ detection method known as PRAISE.\nAs a result, the authors newly identified the RNA substrates and corresponding\nmodification sites for multiple pseudouridine synthases (PUSs). Furthermore,\nthey demonstrated crosstalk among multiple PUS enzymes in the anticodon\nstem-loop region of tRNAs, whereby these enzymes influence each other's\nmodification activity. These findings represent the first comprehensive,\nquantitative, and high-resolution dataset of Ψ modifications and mark a\nsignificant milestone in the field of RNA modification research. Below are\nseveral points for consideration.\nResponse: We thank this referee for the positive comments and the\nconstructive suggestions! To address specific comments:\n1. “Adding a discussion on the physiological relevance of the observed interplay\namong PUS enzymes would greatly enhance the manuscript. For instance,\ncould such interplay reflect a compensatory mechanism that maintains tRNA\nfunctionality in response to the loss of specific modifications?”\nResponse: Thanks for your valuable suggestion. Modifications within the tRNA\nanticodon are responsible for maintaining translational fidelity and efficiency\n(PMID: 18552770). In particular, the physicochemical properties of Ψ stabilize\ntRNA structure, and its presence in the anticodon stem-loop contributes to the\nregulation of translation (PMID: 10446245; PMID: 29507244). Consequently,\nthe loss of a specific modification site may trigger a coordinated response\namong PUS enzymes, reflecting a compensatory mechanism that preserves\ntRNA functionality. We have added several sentences and references to the\ndiscussion section.\n2. “In PUS1 and PUS7 knockout cells, increased Ψ modification at position 39\nwas observed in many tRNAs that contain the respective target sites for these\nenzymes. Do similar changes also occur in tRNAs that do not carry these target\nsites? Discussing this point would be important for understanding potential\ninteractions among PUS enzymes.”\nResponse: Per your suggestion, we performed a detailed analysis of Ψ level\nchanges in tRNAs upon PUS1 or PUS7 KO. First, tRNAs carrying uridine at\nposition 39 were classified into PUS1-targeted (n= 86) and PUS1-non-targeted\n(n= 20) groups based on our original manuscript. Although 9 out of 20 non-\ntargeted tRNAs were initially excluded by stringent filtering criteria (deletion rate\ndifference > 5% and > 2-fold decrease in Ψ level) and thus defined as PUS1\nnon-targets, they exhibited a notably decreased Ψ level at position 27/28\nfollowing PUS1 KO (Reviewer-only Figure 1). Therefore, we redefined these 9\ntRNAs as PUS1 targets. Overall, 95 PUS1-targeted tRNAs showed a globally\nelevated modification level at position 39 in PUS1 KO cells (Extended Fig. 8d).\nIn contrast, for the 11 non-targeted tRNAs, no visible increase in Ψ level was\ndetected at position 39 within the PUS1 non-target group (Extended Fig. 8e).\nSer−CGA−3\n20\nLeu−TAA−3\nSer−TGA−1\nLeu−TAA−1\nCys−GCA−12\nCys−GCA−10\nArg−TCT−4\nWT PUS1KO\nReviewer-only Figure 1 | Ψ levels of 9 modification sites in WT and PUS1 KO, which can\nbe further refined as PUS1-dependent sites.\nFurthermore, tRNAs carrying uridine at position 39 were categorized into two\ngroups: those that contain PUS7 target sites (n= 13) and those that do not (n=\n44). Consistently, in PUS7 knockout cells, 12 PUS7-targeted tRNAs exhibited\n)%(\nlevel\nΨ\nincreased Ψ39 levels, whereas no substantial changes in Ψ modification were\ndetected at position 39 within the PUS7 non-targeted tRNA group (Extended\nFig. 8f-g). These results indicate that the observed Ψ39 elevation arises from\nan interplay that occurs specifically at PUS-dependent substrate sites, rather\nthan from a global compensatory modification effect. We have added the new\nresults and discussion to the revised manuscript.\n3. “To provide a more comprehensive understanding, it would be valuable to\npresent data on the Ψ modification sites and their levels in snRNAs and\nsnoRNAs in RPUSD1–3 knockout cells.”\nResponse: Thank you for the suggestion. To address this point, we analyzed\nΨ modification sites and their levels in snRNAs and snoRNAs from RPUSD1–\n3 knockout cells. Unfortunately, we did not identify any confident RPUSD1–3-\ndependent Ψ sites in these RNA species (Extended Fig. 2a-b). We further\nexamined the PRAISE datasets generated from additional PUS knockout cell\nlines that were included during revision, namely PUS1L, PUS7L, and PUS10,\nand likewise did not identify any Ψ sites in snRNAs or snoRNAs (Extended Fig.\n2a-b). These observations indicated that the majority of Ψ sites in snRNAs are\ninstalled by RNA-guided modification machineries.\nTo further evaluate this possibility, we compared the Ψ sites that showed no\nsignificant response to DKC1 KD with two established snoRNA-targeting\ndatabases (PMID: 16381836; PMID: 27174936). Among these sites, 12 of the\n16 Ψ sites in snRNAs had already been predicted to be snoRNA-guided\n(Supplementary Dataset 6). For the remaining four sites, we performed an in-\nsilico guide prediction analysis using the complete set of human H/ACA box\nsnoRNAs from snoRNABase (PMID: 16381836). We employed stringent\nfiltering criteria in RNAsnoop to minimize false positives: (1) ≥ 4 base pairs in\nthe 5′ duplex, (2) ≥ 4 base pairs in the 3′ duplex, and (3) a total predicted free\nenergy < –10 kcal/mol. Despite these stringent parameters, we successfully\nidentified plausible guiding H/ACA snoRNAs for those 4 Ψ sites\n(Supplementary Dataset 6).\nWe then expanded this analysis to Ψ sites within snoRNAs. Overall, the in silico\nprediction yielded high-confidence guide–target assignments, and 22 of the 25\nΨ sites in snoRNAs could be matched to putative guiding H/ACA snoRNAs.\nThe complete mapping, including base-pairing interactions and predicted\nbinding energies, is presented in Supplementary Dataset 6. As three\nrepresentative examples, we show the predicted interactions between H/ACA\nsnoRNAs and Ψ sites in U5 snRNA and snoRNAs (Extended Fig. 2d).\nCollectively, these results indicated that Ψ installation in snRNAs and snoRNAs,\nwith the exception of the single stand-alone PUS7-dependent site, was\npredominantly mediated by snoRNA-guided modification pathways. These\nresults have been incorporated into the revised manuscript.\n4. “It is useful to summarize all the data (Ψ position, tRNA species, and\ndependency on each enzyme) in a single table (e.g., an Excel sheet).”\nResponse: We appreciate your suggestion. Supplementary Dataset 1 and\nSupplementary Dataset 2 have now been compiled, each listing cyto-tRNAs\nand mt-tRNAs with defined PUS dependencies. These tables report the Ψ\npositions and corresponding tRNA species in WT cells, the modification levels\nin the relevant knockout lines, and the inferred enzyme dependencies.\nOther comments\n1. “The data indicate that Ψ8 is present only in Gly-GCC, but Reference 14\nclaims that Ψ8 is found in fragments derived from Ala, Cys, and Val as well.\nHow do the authors reconcile this discrepancy?”\nResponse: Considering that Ψ8 has not been detected in mature Ala, Cys, and\nVal tRNAs in both our study and previous relevant research (Song, J. et al., Nat.\nChem. Biol. 2020; Cui, Q. et al., Nat. Cancer. 2021; Xu, H. et al., Nat. Methods.\n2025), thus we speculated that Ψ8 observed in tRNA-derived small fragments\nmight be installed after tRNA cleavage.\n2. “For abbreviations such as CMC, the full name must be given on first\nappearance.”\nResponse: Thanks. We have provided the full name of CMC on its first\nappearance.\n3. “Line 118: Does this refer only to cytosolic tRNAs? Please clarify whether\nmitochondrial tRNAs are included.”\nResponse: Thank you for the suggestion. The statement refers only to\ncytosolic tRNAs. We have added the fraction of mt-tRNAs (0.39%) in the\nrevised Figure 1 and rephrased the sentence to make it clearer.\n4. “Line 159: Describe the specific values (total number of sites, etc.) on each\nregion.”\nResponse: We have provided the number of Ψ sites in each tRNA region in\nthe revised manuscript.\n5. “Line 162-163~. The number of Ψ sites at each position is described, but the\nmeaning of this number is difficult. Because Ψ is induced only in uridine, this\nnumber is significantly affected by the ratio of uridine at each position.\nDescribing the uridine number at each position is helpful.”\nResponse: Thanks for the helpful suggestion. We have described the uridine\nnumber at each position in the revised manuscript.\n6. “Line 171: Please clarify this ‘median level’ is calculated from all tRNAs or\nmodified tRNAs only.”\nResponse: The median level is calculated from modified tRNAs only. We have\nclarified this in the revised manuscript.\n7. “Line 242: Instead of “UGG”, would “USG (S = G or C)” or “USR (R = G or\nA)” be more appropriate?”\nResponse: We have revised “UGG” to “USR (S = G or C, R = G or A)”.\n8. “Line 253: ‘Fig. 2c,d’ should be corrected to ‘Fig. 2b,c’.”\nResponse: We have revised the references from “Fig. 2c,d” to “Fig. 2b,c.”\n9. “Line 481: The sentence appears to be unclear or grammatically incorrect.\nPlease revise for clarity.”\nResponse: We have corrected its description in the revised manuscript.\n10. “Line 482: This line is missing or cut off, making it unreadable.”\nResponse: We have corrected to make it readable.\n11. “Line 488: Please provide the DNA sequences used for PCR.”\nResponse: We have provided it in the revised manuscript.\n12. “Lines 497 and 499: Should “8 bp” be corrected to ‘10 bp’?”\nResponse: We have corrected it to “10 bp” in the revised manuscript.\n13. “Line 540: Please clarify what “10” refers to. 10 reads?”\nResponse: This “10” refers to “10 reads.” We have revised it for clarity.\n14. “Figure 1a: Are the values the average of two biological replicates?”\nResponse: The values in Figure 1a represent the average of two technical\nreplicates. We have added a description in the figure legend to clarify this in the\nrevised manuscript.\n15. “Figure 1b: What are the different sites between different repeats? Show\nthe details.”\nResponse: We analyzed the detailed information of a total of 17 non-\noverlapped Ψ sites in two replicates. We found that 12 out of the 17 sites failed\nto overlap due to low modification levels (3%~5%) (Reviewer-only Figure 2),\nand the remaining 5 non-overlapping Ψ sites showed high P values as a result\nof low read coverage in the other replicate and were therefore filtered out\n(Extended Fig. 1e). We have incorporated the results into the revised\nmanuscript.\nReviewer-only Figure 2 | The Ψ levels of the 12 non-overlapping sites in Rep1 and Rep2,\nwhich were filtered out by the 5% deletion-rate cutoff, indicate that these sites exhibit low\nmodification levels in both replicates.\n16. “Figure 1e: The label “arm and loop” should be corrected to “stem loop” as\ndescribed in the text.”\n)%(\nESIARP\nni\nlevel\nΨ\n10\n5\n5% cutoff\n0\np\n1\np\n2\np\n1\np\nResponse: We have corrected this in the revised manuscript.\n17. “Figure 1e legend: Please provide an explanation of the colored bars.”\nResponse: We have rephrased the figure legend in the revised manuscript.\n18. “Figure 1e, 3e, 5d: The heat maps are very hard to read, not only X axis, Y\naxis and signal, please revise them to a more clear version.”\nResponse: We have revised all the heatmaps to improve their clarity.\n19. “Figure 1f: ‘Y8’ is not shown. It would also be helpful to indicate secondary\nstructure information (e.g., anticodon stem loop) to correspond with Fig. 1e.”\nResponse: We have included Ψ8 in Figure 1f and indicated the secondary\nstructure information to correspond with Figure 1e in the revised manuscript.\n20. “Figure 2: Please define the abbreviations “UT” and “NC” in the figure\nlegend. In addition, ensure consistent labeling—for example, “DKC1-NC cells”\nis used in the legend, whereas “DKC1-NC-UT” appears in the figure.”\nResponse: We have defined the abbreviations and ensure consistent labeling.\n21. “Figure 3e,5d: Overall color of this figure is a light pink. Is that OK?”\nResponse: We have adjusted their display format to improve clarity.\n22. “Figure 5g: Is it correct that positions 27 and 35 of Arg-TCT should actually\nbe 28 and 36? Similarly, should position 28 of Met-CAT be 27?”\nResponse: We have corrected it in the revised manuscript.\n23. “Ex. Fig. 1a: “Alkb” should be corrected to “AlkB”.”\nResponse: We have corrected it in the revised manuscript.\n24. “Ex. Fig. 1b: Should “puro” be replaced with “spacer”?”\nResponse: We indeed incorporated a 3’-puromycin modification to prevent\nadaptor self-ligation. To clarify this, we have rephrased the figure legend.\n25. “Ex. Fig. 1b: Should ‘PNK’ be replaced with ‘CIAP’?”\nResponse: We apology for this mistake. We actually used PNK in our study to\nremove the 3’ phosphate from RNA ends, and this has been corrected in the\nMethod section.\n26. “Ex. Fig. 1e: It would be helpful to include a correlation plot between\nPRAISE and BACS.”\nResponse: We appreciate your suggestion. During revision, we performed a\ncomparative analysis of Ψ sites identified in the PRAISE (n = 705) and BACS\ndatasets (n = 607). Because PRAISE reported consecutive Ψs as continuous\nintervals, we define these BACS sites within the same interval as overlapping\nsites. A total of 510 Ψ sites were overlapped in two datasets, covering 94% of\nthe sites in the BACS dataset (Extended Fig. 1i). We further analyzed these\nuniquely detected sites in BACS. Of these, 19 out of 32 displayed a detectable\nsignal in PRAISE as well, though at lower modification levels (2%~5%). The\nremaining 13 sites not detected in PRAISE were attributed to low sequencing\ndepth (n = 2), high background signal (n = 2), different reference versions (n =\n3), or no deletion signal (n = 6) (Extended Fig. 1i-j). Per your suggestion, we\ncompared the detected Ψ levels between PRAISE and BACS. We found that\nmost overlapping sites are highly modified and show strong concordance in\nquantification levels (Extended Fig. 1k). For moderate or low modified sites, the\nΨ levels exhibited variability, which may result from the intrinsic characteristics\nof the chemical labeling strategies employed in the two approaches (Extended\nFig. 1k). We have incorporated the new results and discussion into the revised\nmanuscript.\n27. “Ex. Fig. 1f: A figure legend is missing.”\nResponse: We have added the figure legend of original Ex. Fig. 1f (revised Ex.\nFig. 1h).\n28. “Ex. Fig. 1f–g: Please correct the label to ‘1g–h’.”\nResponse: We have corrected the label to “1l-m” in now revised Ex. Fig. 1.\n29. “Ex. Fig. 2a and 2c: Clearly indicate whether the data show cytosolic or\nmitochondrial tRNAs. The figure legend should specify “left” or “right”\naccordingly. The same applies to Ex. Fig. 5b and 5c, and Ex. Fig. 6b and 6c.”\nResponse: We have revised the figure legends to indicate the tRNA types\nshown in the figures and to correspond the left and right panels with their\nrespective legends. At the same time, we have revised and clarified all figure\nlegend descriptions throughout the manuscript.\n30. “Ex. Fig. 2d, Extended Fig. 5d, Extended Fig. 6a, and Extended Fig. 7a: In\nthese figures, it is not immediately clear that the targeted gene itself is not lost\nin the knockout condition. It would be helpful to add an explanatory note in the\nmain text.”\nResponse: Thanks for your suggestion. We have added the explanation in the\nmain text: “We note that in PUS1-knockout cells, PUS1 mRNA is not lost\ndespite the loss of protein function, likely because the premature stop mutation\ndoes not fully trigger nonsense-mediated RNA decay.”\n31. “Ex. Fig. 3: The DNA sequence is shown, but it is unclear which region of\nthe gene it represents. Please indicate where the indel occurred in the gene\nstructure. Also, figure legend labels (a–f) do not match the figure panels. Was\nWestern blotting performed?”\nResponse: Thank you for pointing this out! We have now added new panels\nfeaturing a schematic diagram of the gene structures (Extended Fig. 4). These\ndiagrams precisely mark the positions where the indels occurred within the\ntargeted exons. We also revised the figure legends to match the corresponding\nfigure panels.\nAdditionally, as you suggested, we performed Western blot analysis. The\nresults demonstrate a complete loss of the RPUSD1 and RPUSD3 proteins in\nthe respective knockout cell lines. These new data have been included as\nExtended Fig. 4g-h.\nPoint-by-point Response to Reviewer #2\nIn this manuscript entitled “Quantitative analysis of small RNA\npseudouridylation reveals interplay of PUS enzymes in tRNA anticodon stem-\nloop”, Liu et al adapted their original PRAISE method to small RNA Ψ profiling\nby incorporating AlkB demethylation to remove the methylation of several tRNA\nmodifications to improve the readthrough in RT reaction, enabling quantitative\nand site-specific profiling of Ψs in small RNA. The method worked very well for\nquantitative profiling the Ψ sites in small RNA including tRNA, snRNA and\nsnoRNA. By perturbing seven PUS enzymes using stable knockdown or\nknockout cell lines, the authors made interesting discoveries including: (1)\npseudouridylation of snoRNAs can be mediated other stand-alone enzyme\nPUS7 in addition to DKC1; (2) quantification of the Ψ fractions in KO and WT\nsamples not only allow the assignment of the tRNA substrates of individual PUS\nenzyme, but also reveal unexpected interplay among PUS enzymes in shaping\nΨ patterns within the ASL; and (3) for the first time the tRNA targets of the\nhuman RluA family were determined the human tRNA targets which differs from\nthat of its yeast counterparts. The manuscript was written very clearly and the\nevidence they showed are adequate to support their conclusions. Their finding\nwill definitely benefit more functional studies of Ψ and PUS enzymes. Therefore,\nthis reviewer supports its publication in Nature Communications after minor\nrevision.\nResponse: We thank this referee for the positive comments and constructive\nsuggestions! To address specific comments:\nThe following are some minor comments:\n1. “Line 107-108, ‘However, specific modifications, including m1A and m3C\nmodifications, can block reverse transcription’. Besides m1A and m3C, how\nabout other modifications in tRNA such as m1G and m22G? Can these\nmodifications also block RT reactions? Can these RNA modifications be\ndemethylated by AlkB treatment?”\nResponse: Thanks for the suggestion. According to previous studies,\nmethylations on the Watson–Crick face of tRNA bases, including m1A, m3C,\nm1G, and m2 G, can impair cDNA synthesis. Specifically, m1A and m3C can be\n2\nremoved by wild-type AlkB, m1G can be demethylated by the D135S mutant,\nand m2 G can be efficiently converted to m2G by the AlkB D135S/L118V double\n2\nmutant (PMID: 28371071). In our study, we employed wild-type AlkB to\ndemethylate the predominant m1A and m3C modifications. To determine the\nreverse transcription (RT) signature of m1G and m2 G using the Maxima H\nand deletion signals. The results revealed that m1G9/37 primarily produced\nstop signals, whereas m2 G26/27 exhibited both mutation and stop signal\n2\n(Extended Fig. 11a).\nTo explore whether Ψ interacts with other modifications in tRNA, we further\nanalyzed the modification levels of RT-detectable sites across wild-type and six\nPUS knockout cell lines in our original manuscript, as well as three additional\nPUS knockout cell lines (PUS1L, PUS7L, and PUS10) included during the\nrevision. Interestingly, when PUS1, which catalyzes Ψ27/28 formation, was\nknocked out, the adjacent m2 G26/27 modification was also significantly\n2\nreduced (Extended Fig. 11b). Because these sites reside within the same tRNA\nisodecoder, this observation suggests a local crosstalk between adjacent\nmodification sites. We have added the new results and discussions to the\nrevised manuscript.\n2. “Line 118-9: “rRNA accounted for 12.3%”. The small RNA fraction should not\nhave much rRNA contamination, is there some rRNA degradation in the course\nof small RNA purifications?”\nResponse: Thanks for pointing this out! We apologize for confusion regarding\nthe rRNA fractions. In our original manuscript, the rRNA references used for\nmapping included only 5S and 5.8S rRNAs, which are components of the small\nRNA population. The majority of mapped rRNA reads corresponded to 5S rRNA\n(12.18%), with a small fraction mapping to 5.8S rRNA (0.04%) (Fig. 1a).\nIn addition, we reanalyzed the WT PRAISE data by mapping them to the 18S\nand 28S rRNA references and found that 0.1% of the reads mapped to 18S\nrRNA and 0.5% mapped to 28S rRNA. Therefore, the observed 12.3% rRNA\npopulation is not attributable to degradation of 18S or 28S rRNAs during RNA\npurification. We have corrected the corresponding panel in Fig. 1a and updated\nthe description accordingly in the revised manuscript.\n3. “Line 123 “In the cyto-tRNAs, 713 and 714 Ψ sites were detected in two\ntechnical replicates, respectively, with 705 overlapping Ψ sites, thus\ndemonstrating high reproducibility…” it would be nice to describe the quality of\nthese libraries, such as the mapping ratio, deduplication rate, insert sizes, and\nthe criteria to call these sites. In addition, for the non-overlapped Ψ sites, are\ntheir fractions lower in general than the detected sites? What might be the\nreasons why these sites failed to overlap?”\nResponse: In response to the request, we have added additional data quality\nmetrics, including the mapping ratio, deduplication rate, and insert size\ndistribution in two replicates of WT samples (Extended Fig. 1d). We have also\nclarified the criteria used for site identification, defined as a deletion rate\ndifference (treated vs. untreated) exceeding 5% and a P value < 0.0001 in\nMethod section. In addition, we further analyzed the detailed information of a\ntotal of 17 non-overlapped Ψ sites in two replicates (Extended Fig. 1e). We\nfound that 12 out of the 17 Ψ sites failed to overlap due to their low modification\nlevels (3%~5%) (Reviewer-only Figure 2), and the remaining 5 non-overlapping\nΨ sites showed high P values as a result of low read coverage in the other\nreplicate and were therefore filtered out.\nReviewer-only Figure 2 | The Ψ levels of the 12 non-overlapping sites in Rep1 and Rep2,\nwhich were filtered out by the 5% deletion-rate cutoff, indicate that these sites exhibit low\nmodification levels in both replicates.\n4. “Line 174-175 “modification level, reaching a median of 83.4%. Ψ35, located\nat the center of tRNA-Tyr anticodon, showed a high modification level of 86.9%”.\nThe detected modification level is calculated based on the deletion rate and the\ncalibration curves?”\nResponse: The modification levels presented in the manuscript were derived\nfrom the median absolute deletion rates measured at specific sites within\nmodified tRNAs. We have added this clarification to the revised manuscript for\nimproved clarity.\n5. “Line 213 “12 high-confidence DKC1-dependent sites among 28 Ψ sites in\nsnRNA”, for the other 16 Ψ sites in snRNA, any of them is sensitive to the KO\nof the other PUS enzymes?”\nResponse: For the remaining 16 Ψ sites that showed no significant response\nto DKC1 knockdown, we further examined their signals in additional PUS\nknockout datasets to assess whether they might be modified by stand-alone\n)%(\nESIARP\nni\nlevel\nΨ\n10\n5\n5% cutoff\n0\np\n1\np\n2\np\n1\np\nsynthases. This analysis included RPUSD1, RPUSD2, and RPUSD3 knockout\ncells from our original manuscript, as well as three additional PUS knockout\nlines, including PUS1L, PUS7L, and PUS10, generated during the revision.\nHowever, our analysis showed that none of these six enzymes was responsible\nfor modifying these 16 sites (Extended Fig. 2a). Thus, the lack of significant\nchange might be attributed to residual DKC1 activity under knockdown\nconditions.\nTo investigate whether the remaining 16 Ψ sites in snRNAs are mediated\nthrough the snoRNA-guided mechanism, we first compared these sites with two\nestablished snoRNA-targeting databases (PMID: 16381836; PMID: 27174936).\nAmong these 16 sites, 6 were independently annotated in both databases, and\nan additional 6 were documented in at least one of them. Therefore, the majority\nof Ψ sites in snRNAs (24/28) have been previously validated or predicted to be\nsnoRNA-guided and catalyzed by the DKC1 (Supplementary Dataset 6). For\nthe remaining 4 sites not previously documented in either database, we\nperformed an in-silico target prediction analysis using the comprehensive\ncollection of human H/ACA box snoRNAs from the snoRNABase (PMID:\n16381836). We employed stringent filtering criteria in RNAsnoop to minimize\nfalse positives: (1) ≥ 4 base pairs in the 5′ duplex, (2) ≥ 4 base pairs in the 3′\nduplex, and (3) a total predicted free energy < –10 kcal/mol. Despite these\nstringent parameters, we successfully identified plausible guiding H/ACA\nsnoRNAs for all 4 Ψ sites (Supplementary Dataset 6). For example, snoRNA\nACA23 is predicted to base-pair with a Ψ site in U5 snRNA (Extended Fig. 2d).\nCollectively, these results further support that Ψ installation in snRNAs\npredominantly relies on snoRNA-guided mechanisms. These results have been\nincorporated into the revised manuscript.\n6. “Line 411-414, “facilitate characterization of additional PUS enzymes, such\nas PUS1L and PUS7L, enabling comprehensive mapping of tRNA Ψ\nmodifications and advancing our understanding of PUS-mediated Ψ formation.”\nIs there any specific reason why the authors don’t include the KO of PUS1L and\nPUS7L in this study?”\nResponse: We appreciate your question. To provide a more comprehensive\nmap of tRNA Ψ modifications, we generated PUS1L-/- and PUS7L-/- HEK293T\ncells and performed PRAISE sequencing to identify potential Ψ targets in both\ncytosolic and mitochondrial tRNAs. Upon PUS1L knockout, no significant\nchanges in Ψ levels were detected in cytosolic tRNAs, whereas seven Ψ39/40\nsites in mitochondrial tRNAs completely disappeared, demonstrating that\nhuman PUS1L specifically catalyzes Ψ in mitochondrial tRNAs (Extended Fig.\n10b). Thus, PUS1L and PUS3 exhibit compartment-specific Ψ activities:\nPUS1L acts in mitochondria, and PUS3 modifies nucleus-encoded tRNAs. We\nalso performed immunofluorescence staining of PUS1L and PUS3 in HeLa\ncells, and consistently found that PUS1L was mainly localized to mitochondria,\nwhereas PUS3 was abundant in the nucleus (Extended Fig. 10c, d).\nUpon PUS7L knockout, we found 25 Ψ sites exhibited a marked decrease in\nmodification levels (Extended Fig. 10e, f). These PUS7L targets occur\nexclusively in type II tRNAs, characterized by their extended variable loop. The\ntargeted tRNAs include Ser-GCT, Leu-AAG, Leu-CAA, and Leu-TAG, and all\nΨ sites were confined to the variable loop region (Supplementary Dataset 1).\nTogether, these results indicate that PUS7L displays substrate specificity\nentirely distinct from its paralog PUS7, suggesting that it may have specialized\nfunctions that warrant further investigation. The new findings have been\nincorporated into the revised manuscript, and the Discussion section has been\naccordingly updated.\nWhile our work was under revision, we noticed that another relevant study was\npublished online reporting the tRNA Ψ targets of several stand-alone PUS\nenzymes based on the BACS method (PMID: 41136621). Hence, we further\ncompared the tRNA substrates of PUS1L and PUS7L using this orthogonal\nmethod and found highly consistent results. We have now cited this study in\nthe revised manuscript.\nPoint-by-point Response to Reviewer #3\n1. “The authors performed a comprehensive analysis of pseudouridylation in\nhuman tRNAs as well as other small RNAs including snoRNAs and snRNAs.\nThe data demonstrate all of the positions of pseudouridylation, while quantifying\nlevels using their previously established method, named PRAISE. They were\nalso able to assign many of the sites to specific pseudouridine synthases\nthrough analysis of cells with knockdown or knockouts of 7 different\npseudouridine synthases. The innovation derives primarily from the large scale\nof data, and this work will serve as a very useful resource for the community\nstudying RNA modifications in tRNAs, snRNAs, and snoRNAs. The authors do\nhighlight some novel findings, including a snoRNA pseudouridine site not\ncatalyzed by DKC1 but by “stand-alone” synthase PUS7, as well as new\nexamples of apparent interplay between different pseudouridines in the\nanticodon stem/loop, in which loss of one site can lead to either increases or\ndecreases of other sites catalyzed by other enzymes. These examples are\nlimited in comparison to the rest of the data presented, however, they are not\northogonally validated by another method and are not considered in regard to\nhypotheses or models that could be tested, which would have extended their\nsignificance for the field. This work represents a great resource for the field, but\nalone does not present a significant new understanding of tRNA biology, although\nthe data have this potential if extended through further analysis or discussion.”\nResponse: Thanks for the suggestion. We agree with the reviewer that\nincluding mechanistic investigation and orthogonal validation of the PUS\ninterplay would further strengthen and extend the significance of our study.\nConsidering that some PUS enzymes possess catalytic-independent functions,\nfor example, TRUB1 can act as a tRNA chaperone (PMID: 27849601), we\nsought to determine whether the observed interplay among PUS enzymes\ndepends on their catalytic activity. Specifically, we transfected PUS1-KO cells\nwith the nuclear-localized wild-type PUS1 isoform (PUS1-WT) and a\ncatalytically inactive mutant (PUS1-D118A), followed by PRAISE sequencing\nof small RNAs (PMID: 18648068). PUS1-WT successfully restored Ψ28 levels\nin Ile-AAT tRNAs, whereas the catalytic mutant failed to rescue the modification\n(Extended Fig. 7f).\nFurthermore, analysis of the neighboring RPUSD1 substrate Ψ30 in Ile-AAT\ntRNA showed that its elevated modification level in PUS1-KO cells was\nmarkedly reduced upon PUS1-WT rescue, thus demonstrating a successful\nrescue of PUS interplay. The incomplete recovery might be due to transient\ntransfection, under which a portion of pre-existing Ile-AAT tRNAs remains Ψ30-\nmodified at steady state. In contrast, the catalytic mutant could not rescue the\ninterplay to the same extent as PUS1-WT (Extended Fig. 7g-h), which suggests\nthat the catalytic activity of PUS1 is essential for its interplay with additional\nPUSs.\nIn line with our findings, reanalysis of published BACS data (PMID: 39349603)\nrevealed a significant increase in the Ψ30 level of Ile-AAT tRNA in PUS1-KO\ncells, providing orthogonal evidence for the interplay between PUS1 and\nRPUSD1 (Extended Fig. 3f). We have added the corresponding descriptions in\nthe main text of the revised manuscript.\n2. “This work showcases the potential of PRAISE to establish the record of\nwhere pseudouridine is and at what levels quantitatively. A limitation is that\ninformation about other modifications is not recorded. This may not influence\ninterpretation of analysis of snoRNAs and snRNAs, because those RNAs are\nnot known to contain a wide variety of modifications, at high density. But for\ntRNAs, there is a lot of lost information in not considering how loss of\npseudouridine may influence other modifications. For example, several\nexamples are known in which loss of pseudouridine in tRNA at specific sites\ncan affect other modifications (Porat 2025 for review). While it remains\nchallenging to record this information, it is suggested that the authors\nnonetheless acknowledge this limitation of the approach in the text, and discuss\nthe possibility that some of the interplay that they observe is “indirect” (in that it\nis happening through one pseudouridine influencing another modification,\nwhich is in turn influencing a different pseudouridine site). To summarize, due\nto the complexity of the tRNA modification landscape, it is difficult to draw\nconclusions for interplay between modifications when only examining one type\nout of more than two dozen possibilities in tRNA.”\nResponse: We appreciate your suggestion. Previous studies have shown that\nseveral RNA modifications, including m1A, m3C, m1G, m2 G, ms2t6A/ms2i6A,\n2\nacp3U, and inosine, can interfere with cDNA synthesis during reverse\ntranscription (RT), thereby inducing misincorporation or termination signatures\nin tRNA sequencing (PMID: 33581077). Since PRAISE sequencing utilizes the\nAlkB protein to demethylate the majority of RT-stop m1A and m3C, we further\nanalyzed the remaining four RT-readable modification types and found that\nthey can be detected as significant RT stops, mutations, or deletions in our RT\nreactions (Extended Fig. 11a).\nWe then assessed whether potential crosstalk exists between Ψ and other\ntRNA modifications by comparing WT and all PUS KO cells in our original\nmanuscript, together with three additional PUS KO cells (PUS1L, PUS7L, and\nPUS10) included during the revision. In general, we did not observe interplay\nbetween Ψ and m1G, ms2t6A/ms2i6A, and acp3U. But interestingly, when PUS1,\nwhich catalyzes Ψ27/28 formation, was knocked out, the modification level of\n17 m2 G26/27 sites showed a preferentially pronounced decrease (Extended\nsites, possibly mediated by structural coupling or cooperative substrate\nrecognition. We have incorporated this finding into the Discussion section in the\nrevised manuscript.\nWe cannot exclude the possibility that the observed PUS interplays are indirect\nand may occur through other tRNA modifications. Future studies may employ\nnanopore direct RNA sequencing to systematically investigate such\nmodification interplays at the single-molecule level (PMID: 37024678). We have\nadded several sentences to the discussion section.\n3. “The quantitative analysis of pseudouridylation site modification is impressive,\nbut sometimes other aspects of quantification need improvement. For example,\nin the analysis presented in lines 227–236, it describes how 48/73 sites are\nDKC1-dependent. Then after analysis of three mutants, it is concluded that one\nsite is PUS7-dependent. Was this site one of the 25 that were not DKC1-\ndependent? Presumably? But what about the other 24? Are there any ideas for\nwhat those could be catalyzed by? Given that this result was stated in the\nabstract, there needs to be a more quantitative description of the number of\nsites in snoRNAs that are DKC1-dependent, and the number of new sites that\nhave been revealed to be depended on a stand-alone synthetase, which from\nthis study appears to be one (SNORD65). Moreover, because of the apparent\nsignificance of this result—it is mentioned in the abstract—the authors can\nconsider whether they would like to demonstrate it by an orthogonal approach.\nGiven that the entire study basically consists of PRAISE assays, this type of\nresult is of the type that may most benefit from an orthogonal validation.”\nResponse: Thank you for pointing this out! We have conducted additional\nanalyses and experiments to address the concerns you raised, as detailed\nbelow:\n1) Orthogonal validation of the stand-alone Ψ site in SNORD65. Reanalysis of\nthe BACS dataset confirmed that this site is PUS7-dependent, with the Ψ level\nat position 47 decreasing from 51.6% to 0.56% in PUS7 knockout cells\n(Extended Fig. 2c).\n2) Exploration of the remaining 25 DKC1-independent sites (one site was\nmissed in our original manuscript). We analyzed Ψ level in snoRNA from PUS1,\nPUS7, TRUB1, RPUSD1, RPUSD2, and RPUSD3 knockout cells, and also\nincluded newly generated knockout cell lines for PUS1L, PUS7L, and PUS10.\nHowever, the results showed that none of the stand-alone PUSs are\nresponsible for modifying these 24 sites (Extended Fig. 2a, b). Hence, we\nspeculated that these sites may also be installed through RNA-guided\nmechanisms, as residual DKC1 activity in the knockdown cells may have been\nsufficient to maintain basal Ψ levels.\n3) Prediction of guiding H/ACA snoRNAs for the remaining 25 DKC1-\nindependent sites was performed. We performed an in-silico target prediction\nanalysis based on the snoRNA database (PMID: 16381836). For each\npseudouridine site, a 21-nt window containing 10 nt upstream and downstream\nof the modified position was extracted to generate the target set. The 3′\ncatalytic hairpin regions of all H/ACA box snoRNAs were obtained using\nRNAfold from the ViennaRNA Package (PMID: 18424795) and used as the\nsnoRNA query library. This analysis yielded high-confidence guide–target\nassignments, and 22 of the 25 Ψ sites in snoRNAs were matched to potential\nguiding H/ACA box snoRNAs (Supplementary Dataset 6). A representative\nexample showing the predicted interaction between ACA20 and the\npseudouridine site at position 46 in SNORD21 is provided in Extended Figure\n2d. The manuscript, including abstract and main text, has been revised\naccordingly.\n4) Dynamic changes in the Ψ level of snRNA and snoRNA under stress\nconditions. Since previous studies in yeast have reported stress-induced\npseudouridylation of U2 snRNA and mRNA (PMID: 21131909, PMID: 25219674),\nwe sought to determine whether similar dynamic changes occur in\nspliceosomal snRNAs or snoRNAs in human cells. To address this, we\nsubjected HEK293T cells to heat shock and oxidative stress (H₂O₂ treatment)\nand performed PRAISE to quantify Ψ levels. In snRNA, we found that the level\nof Ψ89 in U2 snRNA was significantly increased following heat shock, whereas\nthe levels of other Ψ sites in snRNA remained largely unchanged (Extended\nFig. 2e). We did not observe significant changes in snRNA pseudouridylation\nunder H₂O₂ treatment (Extended Fig. 2e). In contrast, snoRNA exhibited more\npronounced stress-responsive dynamics. We identified four sites that showed\nsignificant decreases under heat shock and two sites that decreased under\nH₂O₂ treatment (Extended Fig. 2f). Notably, Ψ95 in SNORA61 sites consistently\ndecreased under both stress conditions, suggesting that this site may represent\na conserved stress-responsive Ψ site in human snoRNAs.\nAll of these new analyses and results have now been incorporated into the\nrevised manuscript. We sincerely thank you for your valuable suggestions,\nwhich have significantly improved the rigor and completeness of our work.\n4. “The description of changes on mt-tRNAs resulting from the absence of a\nTRUB1 were interesting but would benefit from the addition of a description of\nthe possible role of TRUB2. It would also be beneficial to know whether the\ncytosolic tRNAs that were not changed by TRUB1 deletion lacked the motif\ndescribed in line 275. Or differed in other features. The fact that a relatively\nsmall number of these isoacceptors were affected, and to incomplete levels is\ninteresting, and more discussion of PUS10 vs. TRUB1 would improve impact\nof understanding these interesting these results.”\nResponse: Thank you for the suggestion. To address the concerns, we\nperform additional experiments and analyses as detailed below:\n1) During revision, we performed immunofluorescence staining of TRUB2 in\nHeLa cells and found that TRUB2 is predominantly localized to mitochondria\n(Extended Fig. 10a). Previous studies have shown that TRUB2 is essential for\nmitochondrial protein synthesis (PMID: 27974379). However, in our original\nmanuscript, we assigned all Ψ55 sites in mt-tRNA to TRUB1 and all Ψ sites in\nmt-mRNA to PUS1 (PMID: 36997645), thereby excluding a potential\ncontribution of TRUB2 based on our quantitative analyses. Hence, we have\nadded the sentences to the revised manuscript discussing the possibility that\nTRUB2 may exert a non-catalytic role in mitochondria.\n2) We performed motif analysis for Ψ55 on the cytosolic tRNAs that did not\nshow significant changes upon TRUB1 deletion (i.e., the deletion rate did not\ndrop to less than half of the original level) in our original manuscript. However,\nwe found that these sites displayed a strong “GWUCR (W = U or A; R = A or\nG)” motif (Reviewer-only Figure 3a), which is highly similar to the “GUUCR”\nmotif observed at TRUB1-dependent Ψ sites (Reviewer-only Figure 3b).\n2.0\n1.5\n1.0\n0.5\n0.0\n1 2 3 4 5\nPosition\nReviewer-only Figure 3 | a. Sequence motif TRUB1-independent Ψ55 sites. b. Sequence\nmotif TRUB1-dependent Ψ55 sites.\n3) To investigate the potential redundancy between TRUB1 and PUS10 in\ncytosolic tRNAs, we generated a PUS10-/- HEK293T cell line and performed\nPRAISE analysis on small RNAs. Our results showed that Ψ54 sites in tRNA\nwere completely lost upon PUS10 knockout. In addition, 21 sites in Ψ55\nexhibited substantial decreases in PUS10-/- cells, including those in Ala-AGC,\nVal-AAC, Val-TAC, and Val-CAC tRNAs (Extended Fig. 9e). 14 sites in Ψ55\ndisplayed modest decreases and were affected by both PUS10 and TRUB1\nknockout (Extended Fig. 9f). Together, these results revealed functional\nredundancy between PUS10 and TRUB1 at cytosolic tRNA Ψ55 sites.\nWe have incorporated these newly obtained results, along with a revised\ndiscussion of the roles of TRUB1 and PUS10, into the manuscript to improve\nthe clarity and interpretation of our findings.\n5. “In the case where there is “interplay” between different modified sites in the\nanticodon loop, how is “labeling efficiency” of a given site affected when there\nis a strong decrease in a neighboring site? Presentation of data that would\nrefute such an artifactual model for changes would be appreciated.”\nResponse: We thank the reviewer for raising this concern. To address it, we\nchemically synthesized two ~70-nt spike-in RNA oligos: one containing two Ψ\nsites within a sequence resembling the anticodon stem-loop (ASL) of Ile-AAT\nstiB\nTRUB1-independent sites motif\n2.0\n1.5\n1.0\n0.5\n0.0\n1 2 3 4 5\nPosition\nstiB\na b\nTRUB1-dependent sites motif\ntRNA, and, for comparison, another containing a single Ψ site (Extended Fig.\n3d). The fully modified RNAs were mixed with their corresponding unmodified\ncounterparts in ratios ranging from 0% to 100%, and the deletion signals were\nquantified by PRAISE. The expected and experimentally determined Ψ levels\nexhibited excellent quantitative agreement for the neighboring Ψ33 and Ψ37\nsites in oligo-ΨΨ (R2 = 0.99). Specifically, the spike-in with two 100% Ψ sites\nshowed 86% and 82% modification levels after PRAISE detection, whereas the\nsingle-Ψ spike-in showed 84%, indicating that adjacent Ψ sites do not\ncompromise labeling efficiency. Hence, the results demonstrate that the\ninterplay observed by PRAISE is not due to artificial factors but represents a\ngenuine biological phenomenon. We have incorporated the above data and\nanalysis into the revised manuscript (Extended Fig. 3e). Thanks for the helpful\nsuggestions!\n6. “Pathogenic variants of PUS1 and PUS7 are included in the introduction, but\nnot in the discussion. Given that this study found several examples where PUS1\nor PUS7 catalyzed pseudouridine sites impact other sites, it would be\nreasonable to speculate how these interactions may be involved in a disease\nmechanism. In the case of PUS1 disease variants, several are known to have\na reduction/loss of catalytic activity, which could possibly cause altered\ninteractions with other pseudouridine sites compared to a PUS1 knockout.\nIn summary, connecting the findings of this study to pathogenic variants\nwould emphasize the broader importance of the results.”\nResponse: We appreciate the insightful suggestions from the reviewer. In\nresponse to this comment, we investigated how pathogenic variants of PUS1\nmay influence the observed PUS interplay. Based on a reported disease-\nassociated nonsense mutation E220X (PMID: 21686963), we introduced the\ncorresponding nonsense variant E192X into the nuclear-localized PUS1\nisoform 2 and performed rescue experiments in PUS1 knockout cells. Wild-type\nPUS1 and a catalytic-dead D118A mutant were included as controls. Structural\nmodeling using AlphaFold3 showed that the premature termination variant lost\nboth the RNA-binding and catalytic domains, whereas the catalytic-dead\nmutant retained an intact RNA-binding pocket (Extended Fig. 7a-c).\nAt the RNA level, the nonsense variant displayed expression comparable to\nwild-type PUS1 and the catalytic-dead PUS1 mutant. However, western blot\nanalysis detected only trace amounts of the E192X protein, suggesting that this\ntruncation variant triggers rapid proteasomal degradation (Extended Fig. 7e).\nPRAISE sequencing of small RNAs further revealed that wild-type PUS1\nefficiently restored the modification level of Ψ30 in Ile-AAT tRNA. Neither the\ncatalytic-dead mutant nor the nonsense variant showed measurable restoration\nof the interplay phenotype compared with wild-type PUS1. These results\ndemonstrate that pathogenic variants lacking catalytic activity behave similarly\nto PUS1 knockout. All corresponding analyses and results have been included\nin the revised manuscript (Extended Fig. 7f-h).\nSince the missense mutation R116W is located within the conserved catalytic\ncenter and has been predicted to reduce catalytic activity (PMID: 15108122,\nPMID: 28045575), it may influence PUS interplay in a manner distinct from\nPUS1 knockout. We have added a discussion on the potential differential\nimpacts of human pathogenic PUS1 variants on PUS interplay in the revised\nmanuscript.\n7. “The paper title and results section are centered on tRNA. The section on\nsnoRNA and snRNA, while adding value to the study, feels out of place. It is\nsuggested to move this section to the end to maintain the focus on tRNA that\nappears to be the intent of the authors. And to reframe the study as a focus on\ntRNA, but that additional analyses of these other classes of RNA added further\nnovelty to the study.”\nResponse: We appreciate the reviewer’s interest in this part of the study. As\npseudouridylation in snRNAs and snoRNAs was also highlighted by the other\nreviewers, we expanded this section accordingly. To address these comments,\nwe substantially expanded our analyses by validating the stand-alone Ψ site in\nSNORD65 using the BACS dataset, examining the remaining DKC1-\nindependent sites across nine PUS knockout cell lines (including newly\ngenerated PUS1L, PUS7L, and PUS10 knockouts), and performing\ncomprehensive in silico predictions to identify their potential guiding H/ACA\nsnoRNAs. In addition, we assessed stress-responsive pseudouridylation by\nanalyzing snRNA and snoRNA Ψ levels under heat shock and oxidative stress\nconditions. We believe that these data have substantially strengthened this part\nof the study. All generated results and analyses have been incorporated into a\nnew figure in the revised manuscript (Extended Fig. 2).\n8. “The analysis of mRNA and tRNA expression are not sufficiently explained.\nWhy were these done? How do they facilitate interpretation of the PRAISE data?\nSome examples from the literature should be provided that defend the idea that\nloss of one PUS enzyme or other RNA modification enzyme has been shown\nto lead to changes in another enzyme, and to a significant enough degree that\nit would affect the activity of that second enzyme. Otherwise, this analysis\nseems to yield limited value for the space that it occupies, and without further\ndiscussion to emphasize its importance.”\nResponse: Thanks for raising this important question. In our original\nmanuscript, we analyzed the expression levels of other PUS enzymes upon\nknockout of a single PUS gene. This was motivated by previous studies\nrevealing crosstalk between RNA m6A modification and histone modifications,\nwhere m6A can have either positive or negative effects on the expression of\nhistone methyltransferases. Upon Mettl3 knockdown, the expression levels of\nhistone-modifying enzymes changed, leading to corresponding alterations in\nhistone modification levels and resulting in dysfunctional neuronal and\nepidermal development (PMID: 31154015; PMID: 37656787). To further\nhighlight the importance of analyzing the expression of related enzymes, we\nhave added corresponding descriptions and new references in the revised\nmanuscript.\n9. “A major focus of this study is the interplay between pseudouridine sites in\nthe anticodon stem-loop, of which they found several interesting examples.\nHowever, there are many tRNAs that are modified by PUS1 and do not have\nan increase in Ψ levels at U30/31 (Fig. 3e). Are there any shared features (i.e.\ntRNA sequence, structure, etc.) that may influence whether an interaction\nbetween Ψ27/28 and Ψ30/31 is present? In lines 386-387, the authors state\nthat interactions between PUSs primarily occurs in the anticodon stem-loop,\nwhich is somewhat contradictory to the data presented in Fig. 3e and 5d that\nshow several examples of interactions between Ψ sites in other regions of the\ntRNA. The mechanisms proposed in lines 429-439 could benefit from additional\nreferences to structural/biochemical literature on PUSs.”\nResponse: We thank the reviewer for the constructive and insightful comments.\nWe have carefully addressed the concerns as follows:\n1) To address the first point, we categorized tRNAs carrying uridine at positions\n30 or 31 into PUS1-interplay and non-PUS1-interplay groups and performed\nmotif analysis for each group. Those showing increased Ψ levels were enriched\nfor a “RGUSY (R = G or A, S= G or C, Y= C or U)” motif, whereas tRNAs without\nincreased Ψ levels exhibited a “CKUMK (K = G or U, M = A or C)” motif\n(Extended Fig. 3h), indicating that distinct sequence features might account for\ntheir differential susceptibility to PUS1-mediated interplay.\n2) Description of the location of PUS interplay. While our data show that the\nmost prominent examples of PUS interplay occur within the anticodon stem-\nloop (ASL), additional instances can also be observed in other regions,\nparticularly as illustrated in Figure 3e. Accordingly, we have revised the\ndiscussion and added a statement noting that potential interplay in other tRNA\nregions warrants further investigation. We greatly appreciate this suggestion,\nand we believe this will broaden the significance of our study.\n3) We have added additional references on the structure of PUS1 and TRUB1\nenzymes (PMID: 23707380, PMID: 11779468).\n10. “Some of the figure legends are very short and do not provide the complete\ninformation necessary to interpret the plot. For example, there is insufficient\nexplanation of the heat maps shown in 1e, 3e and 5d/f. For example, it is not\nclear the order that the isodecoders are labeled, and in general there just needs\nto be more detail in how these are presented.”\nResponse: Thanks for pointing this out! We have reorganized the panels in Fig.\n1e, Fig. 3e, and Fig. 5d/f, and added the iso-acceptor information on the Y-axis\nto improve clarity. The iso-decoders are arranged in ascending numerical order.\nWe have revised the corresponding figure legends.\nMinor comments:\n1. “Figure 3F, 4J: the labels caused some confusion. Is this a standard\nnomenclature for listing a tRNA isodecoder and its corresponding position?\nSimply adding the word “position” prior to the last number would distinguish the\ngene name and the position in the tRNA of those mapped reads. Alternatively,\nthe data can be split by position, for example using brackets on the right side\nof the plot to group all “position 27” rows together.”\nResponse: We have split the tRNA isodecoder and its corresponding position,\nand then put all “position” rows to the right side of the plot in the revised\nmanuscript.\n2. “line 427–429: these citations are incomplete in presenting what is known\nabout how psi55 influences m1A58. There are numerous cases known in yeast\nin which the influence is positive and not negative (Barraud et al 2019, Shaw\nand Thomas et al 2024, reviewed in Porat 2025).”\nResponse: Thanks for this suggestion. We have now cited this comprehensive\nreview (PMID: 40044218) in the revised manuscript.\n3. “lines 251–252: the meaning of this statement is not clear. Why would a\ndifference in tRNA expression between two samples cause difficulties in\nassignment of the modified positions in tRNA?”\nResponse: We have corrected this statement in the revised manuscript to\nmake it clearer.\n4. “Figure 5D: The figure legend describes the plot as a histogram, when 5D is\na heatmap. Descriptions for 5D and 5E may have been swapped in legend.”\nResponse: Thanks for pointing out this mistake. We have corrected the legend\nof Figure 5D and 5E.\n5. “Figure 5G: this illustration and legend need to be revised for clarification.\nThere are no arrows being shown, and the identity of starting and terminal\nnodes is not clear and thus direction of change. The first two simply have solid\nlines connecting the two sites.”\nResponse: We have corrected the presentation of the illustration and legend\nin the revised manuscript to make it clearer.\n6. “The Github link does not work and the GEO accession number states that\nthe data will not be public until 2027. Therefore, it is not possible to reproduce\nthe analyses with these current circumstances.”\nResponse: We have uploaded the bioinformatics analysis pipeline to GitHub\n(https://github.com/LiuWenQing657/PsiUInterplayAnalysis) and ensured that it\nfunctions properly. In addition, we have made the sequencing data publicly\naccessible through GEO\n(https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE299274).\nAdditionally, we have uploaded the newly generated raw sequencing data\nproduced during the revision process to the National Genomics Data Center\n(NGDC) BioProject and made them publicly available for download\n(https://ngdc.cncb.ac.cn/bioproject/browse/PRJCA050602).\nREVIEWERS' COMMENTS\nReviewer #1 (Remarks to the Author):\nThis revised version has adequately addressed this reviewer's comments.\nFor the discordant sites between different replicates, the authors performed a\ndetailed analysis. These sites all either have low modification levels or show\nvariability between replicates due to low coverage. I think these are tolerable sources\nof noise and they do support the robustness of PRAISE.\nIn addition, the authors have now added an overlap analysis with BACS, and most of\nthe sites are shared between the two methods. This again supports the reliability of\nPRAISE.\nThis reviewer only has one comment. No. 22 has not really been addressed. In Fig.\n5e, the interplay in Arg-TCT is shown between Ψ28 and Ψ36, but in Fig. 5g it\nsuddenly becomes an interplay between Ψ27 and Ψ35 in Arg-TCT. Likewise, in the\nheat map of Fig. 3e, the interplay in Met-CAT is clearly shown between Ψ27 and\nΨ31, but in Fig. 5g it is labelled as Ψ28 and Ψ31.\nWith this minor revision, this reviewer fully supports the publication of the manuscript.\nResponse: Thank you for your correction. In Fig. 5g, we have revised the interplay\nsites on Arg-TCT to Ψ28 and Ψ36, and the interplay sites on Met-CAT to Ψ27 and\nΨ31. We greatly appreciate your support for our publication in Nature\nCommunications.\nReviewer #2 (Remarks to the Author):\nThe authors have addressed all my questions and I support its publication in Nature\nCommunications.\nResponse: Thank you for your thorough review and positive feedback. We greatly\nappreciate your support for our publication in Nature Communications.\nReviewer #3 (Remarks to the Author):\n(Since this is the second round of review, below I only provide comments in relation\nto the requested/suggested changes from the first round, and the author response\nand manuscript revision.)\nMajor comments:\nThe authors thoroughly addressed many of the reviewer comments by performing\nadditional experiments and further analysis of existing data. Therefore, this reviewer\nsupports the publication of this manuscript in Nature Communications pending\nconsideration of a few recommendations below.\nThe authors compare their overall number of Ψ sites to a similar study that was\npublished during revision (PMID: 41136621). It would be beneficial for the authors to\nsummarize (in a few sentences) the major similarities and differences between these\nstudies in the Discussion. The authors have already done an admirable job of\ncomparing the datasets given the timing of when the other study became available,\nbut a few additional sentences in the Discussion will be helpful to readers to see the\noverall big picture. It is beneficial for the community that these two studies are now\navailable and it emphasizes the importance of this study.\nResponse: Thank you for your constructive suggestion. We have now added a\nsummary of the major similarities and differences between our study and this study\n(PMID: 41136621) in the Discussion. We believe this revision strengthens the\nmanuscript as suggested.\nMinor comments:\nLine 194: space within word “c onserved\"\nResponse: Thank you for your correction. We have deleted the space within word “c\nonserved”.\nLines 286-287: “These results demonstrated that human PUS1 was a multifunctional\nenzyme that targets both cyto- and mt-tRNAs.” The mitochondrial isoform of PUS1\nwas previously established to modify mt-tRNAs (PMID: 15772074, PMID: 21686963).\nThe data presented in this manuscript provides a comprehensive list of mt-\ntRNAs/sites modified by PUS1 but the language can be adjusted to better reflect\nprior discoveries.\nResponse: Thank you. We have revised the description in the manuscript per your\nsuggestion.\nLines 389-391: “we performed rescue experiments in PUS1 knockout cells using\nwild-type PUS1 together with the catalytic-dead D118A mutant and the disease-\nassociated nonsense variant E192X.” Changing “together with” to “or” would clarify\nthis statement. A brief description of how the rescue experiments were performed\nwould be helpful (i.e. “WT PUS1, (etc), constructs were transfected into PUS1 KO\ncells..”). The specific PUS1 isoform used in these experiments can also be included\nhere, since the position of the mutations depends on this detail.\nResponse: Thank you. We have revised the description as your suggestion. A brief\ninstruction of the experimental procedure and the specific PUS1 isoform used have\nnow been included in the Results section. Additional methodological details have also\nbeen provided in the Methods section.\nLine 394: “In contrast, neither D118A nor E192X rescued either modification.” Since\nΨ30 levels are expected to decrease if PUS1 catalytic activity was restored, it would\nbe more clear to specify that neither mutant “rescued Ψ28 modification or the\ninteraction with Ψ30”\nResponse: Thank you. We have revised the description in the manuscript as your\nsuggestion.\nExtended Data Fig 1k: Description of the blue and grey dashed lines are not included\nin the figure legend. What is the p-value for the reported correlation?\nResponse: Thank you for your suggestion. We have added these information in the\nlegend and the plot of Extended Data Fig 1k.\nExtended Data Fig 11a: Are these boxplots showing wild-type data? Add labels.\nResponse: Thank you for your suggestion. We have added the labels of these\nboxplots in Extended Data Fig 11a.\nExtended Data Fig 11b: There is no label on the y-axis.\nResponse: Thank you for your suggestion. We have added the label of y-axis in\nExtended Data Fig 11b.\nReviewer #3 (Remarks on code availability):\nI have not reviewed it but was able to access it following author instructions. Thank\nyou to the authors for making it public.\nResponse: Thank you for your interest in our work and for taking the time to access\nand review our codes. We are glad that the publicly available materials were helpful.\nReviewer #4 (Remarks to the Author):\nI co-reviewed this manuscript with one of the reviewers who provided the listed\nreports. This is part of the Nature Communications initiative to facilitate training in\npeer review and to provide appropriate recognition for Early Career Researchers who\nco-review manuscripts.\nResponse: Thank you for participating in the co-review process. We appreciate your\ncollaborative effort in evaluating our manuscript and providing valuable feedback."}]

Decision letter