IntroductionThe 5′-terminus of eukaryotic mRNA is characterized by a prevalent cap structure composed of 5′,5′-triphosphate-linked 7-methylguanosine (m7G). This canonical capping exerts profound impacts on modulating various functions of mRNA species, such as stability, splicing, polyadenylation, and translation1,2. Recent studies have revealed that nucleotide-containing metabolites, including nicotinamide adenine dinucleotide (NAD), flavin adenine dinucleotide (FAD), 3′-dephospho-coenzyme A (dpCoA), uridine diphosphate glucose (UDP-Glc), uridine diphosphate N-acetyl glucosamine (UDP-GlcNAc), and dinucleoside polyphosphates, can be incorporated into the RNA 5′-terminus as noncanonical initiating nucleotides (NCIN), resulting in NCIN-RNA3,4,5,6,7. Since the first discovery of NAD-capped RNA in bacteria8, NAD-RNAs have been reported in all kingdoms of life9,10,11,12,13,14,15,16. Although studies have indicated the presence of FAD, dpCoA, UDP-Glc/GlcNAc, dinucleoside polyphosphate caps across different phyla13,17,18,19,20, comprehensive profiling of NCIN-RNAs at the epitranscriptome level is only beginning to be revealed.The current profiling strategies are limited to detecting specific metabolite caps9,18,21,22,23,24,25,26. NAD captureSeq9, along with our recently developed methods, ONE-seq24 and DO-seq26, employ chemo-enzymatic reactions to label NAD-capped RNAs with biotin, catalyzed by adenosine diphosphate-ribosyl cyclase. These biotinylated RNAs are then captured by streptavidin beads and subjected to high-throughput sequencing to identify NAD-RNAs. Yet, these methods only reflect the relative enrichment for a specific type of metabolite-capped transcript and cannot quantify the ratio between metabolite-capped and m7G-capped forms. Decapping-based strategies, such as CapZyme-seq27, utilize processing enzymes like NudC to target NCIN-capped RNAs, yielding RNAs with a 5′-monophosphate. Subsequently, barcoded RNA adapters are ligated to the 5′-end of these RNAs, followed by adapter-specific PCR enrichment and high-throughput sequencing. However, a notable limitation of applying CapZyme-seq in eukaryotes arises from contaminating signals of m7G, the dominant cap type28. Furthermore, currently available methods lack the ability, including experimental and computational procedures, for the epitranscriptome-wide percentage assessment of NCIN capping at the transcript resolution.Aging is characterized by physiological declines accompanied with metabolic and transcriptional alterations29,30, but how these processes are integrated into the regulation of aging remains elusive. Studies in rodents and humans have shown that the levels of NAD and FAD decline with age in multiple tissues, including the liver31,32,33. Additionally, our recent findings reveal prominent age-associated alterations in NAD capping in mouse liver and human peripheral blood24,34, highlighting the potential role of metabolite-modified RNAs during aging. Cellular metabolism and gene expression are fundamental biological processes, making it crucial to explore how metabolite-modified RNAs, which coordinate metabolic and transcriptional regulation, are modulated during the adult life cycle and consequently, their physiological impact on the progression of aging. However, the absence of quantitative methods to determine epitranscriptome-wide NCIN caps with stoichiometric information impedes investigations into the contribution of NCIN capping to a broad spectrum of biological processes.In the present study, we demonstrate that the CompasSeq analytical platform, combined with experimental and computational procedures, enables a comprehensive and quantitative assessment of NCIN capping levels in eukaryotes. We further devise an orthogonal and independent approach, the quantitative exoribonuclease reduction assay, to validate newly identified NCIN-RNAs and their capping ratios. Utilizing CompasSeq, we uncover previously unknown features of metabolite-modified RNAs in adult and aging mouse livers.ResultsThe workflow of CompasSeqWe designed the CompasSeq method for Comprehensive percentage assessment of NCIN-RNAs through enzymatic decapping followed by 5′-end adapter ligation and high-throughput Sequencing, an analytical platform that integrated both experimental and computational frameworks (Fig. 1).Fig. 1: The workflow of CompasSeq.a CompasSeq implements a parallel experimental pipeline to simultaneously capture NCIN-capped RNAs and all-capped RNAs (NCIN- and m7G-capped RNAs). Meanwhile, CompasSeq implements a dual spike-in strategy to quantify NCIN-capped RNAs. A set of synthetic NCIN-RNAs with known capping ratios, is introduced at the beginning to serve as internal controls. Adapter-pre-ligated Drosophila RNA, is added at the step of cDNA library construction for subsequent computational analysis. b Diagram illustrating the workflow of CompassAnalyzer. CompassAnalyzer first applies stepwise normalizations to recover the intrinsic RNA abundance distribution between NCIN-capped and all-capped RNAs from sequencing data. It then employs empirical Bayes shrinkage to enable quantification of NCIN capping ratios at the transcript level.Full size imageFirst, we subjected total RNA transcripts to cap-specific enzymes. To selectively capture NCIN-capped RNAs, we began with yDcpS, a decapping enzyme that hydrolyzes m7G-capped mRNA, followed by CIAP treatment to eliminate 5′-phosphate residues. We then introduced the bacterial RNA-decapping enzyme NudC, which hydrolyzes the diphosphate group in the cap structures of NAD, NADH, FAD, and dpCoA, ensuring that only NCIN-RNA yields 5′-monophosphate ready for ligation. To capture both NCIN- and m7G-capped RNAs, we initially treated total RNA with CIAP to remove pre-existing 5′-phosphate moieties. We then applied RppH, a pyrophosphatase, to hydrolyze pyrophosphate from the 5′ end of triphosphorylated RNA, together with NudC, to generate 5′-monophosphate RNAs from all-capped transcript forms. Second, we performed ligation between 5′ oligoribonucleotide adapters and 5′-monophosphate RNAs derived from either NCIN- or all-capped transcripts, followed by polyA-selected reverse transcription and high-throughput sequencing. Of note, we deliberately introduced an exogenous RNA spike-in as an internal control essential for subsequent computational analysis (Fig. 1a). Third, we proceeded with a computational framework composed of stepwise normalizations and empirical Bayes shrinkage. Stepwise normalizations, including global scaling using exogenous RNA spike-ins, linear adjustment, and the RUV strategy35, were applied to recover the intrinsic RNA abundance distribution between NCIN-capped and all-capped RNAs from sequencing data. Empirical Bayes shrinkage was then employed to enable quantification of NCIN capping ratios at transcript resolution (Fig. 1b).Development of CompasSeq experimental proceduresWe generated 38 nt RNA spike-ins via in vitro transcription, using adenine exclusively as the initial nucleotide in the DNA template to synthesize NAD, dpCoA, and FAD-capped RNAs (Fig. 2a). As anticipated, capped RNAs (NAD and dpCoA) exhibited higher molecular weight than uncapped forms (ppp-RNA) as observed by UREA-PAGE gel electrophoresis (Fig. 2a). Consistent with previous reports22,36, two RNA species for FAD-RNA were observed, with upper band being intact FAD-capped RNA and the lower band resulting from the loss of the flavin moiety together with one phosphate group (Fig. 2a). To corroborate the quality of the RNA spike-ins, we employed boronate affinity electrophoresis. Boronyl groups from boric acid form relatively stable complexes with vicinal diols, naturally occurring in the m7G of the m7G cap and the nicotinamide riboside of the NAD cap, leading to gel retardation during electrophoresis37. Indeed, only RNA with m7G and NAD caps, but not ppp-RNA, exhibited electrophoretic mobilities, indicating the presence of these specific caps in the RNA spike-ins (Fig. 2b).Fig. 2: Development of CompasSeq.a In vitro transcription with ATP, m7GpppA, NAD, FAD, or dpCoA as the initiating nucleotide. Sequence of the DNA template with the T7 promoter (Top). The “A” in red was the transcription start site and the coding sequence has been underlined. In vitro-transcribed RNAs with different 5′ end caps were resolved in an 8% PAGE gel (Bottom). The arrows indicated the intact FAD-capped species. b Boronate affinity electrophoresis showing in vitro-transcribed 38 nt RNA, including ppp-, m7G-, and NAD-capped RNAs. c UREA-PAGE gels showing products of processing of NAD- and m7G-RNA 5′-ends by NudC, dRai1, RppH, and yDcpS. d UREA-PAGE gel showing products of processing of FAD- and dpCoA-RNA 5′-ends by NudC. e Analysis of the effect of RNA secondary structure on yDcpS decapping. All m7G-capped RNA hairpins and linear 30-mers (1 µg each) were treated by 9 µM yDcpS for 1 h and then resolved in a 12% PAGE gel. f After NudC treatment, NCIN-RNA, such as NAD-RNA (38 nt), could be ligated with the 5′-adapter (27 nt), as shown by the shift of an upper band in the UREA-PAGE gel (left). m7G-capped RNA (38 nt) after being processed by yDcpS and CIAP, could not be ligated with the 5′-adapter (middle). RppH-processed m7G-RNA, on the other hand, contained a 5′-phosphate group for adapter ligation (right). Source data are provided as a Source Data file.Full size imageAssisted by 38 nt RNA spike-ins, we examined the specificity of decapping enzymes (Fig. 2c). NudC, known for its ability to hydrolyze diphosphate but not triphosphate linkages, effectively removed the NAD cap but also exhibited non-specific cleavage of the m7G cap upon extended treatment (Fig. 2c, d). We also evaluated exoribonucleases from the DXO/Rai1 family, known for their decapping activity27,38. We purified the dRai1 protein, a Drosophila homolog of Rai1, and found that it selectively degraded NAD-capped RNA while leaving m7G-capped RNA intact (Fig. 2c). However, dRai1 also caused additional nucleotides removal due to its exoribonuclease activity (Fig. 2c). In contrast, yDcpS specifically targeted m7G cap, while RppH removed both m7G and NAD caps, lacking specificity (Fig. 2c). It has been reported that RNA secondary structure would hinder the activity of yDcpS39. We generated m7G-capped RNAs with polymorphic hairpins and 3′ or 5′ overhangs to test whether excess yDcpS (~9 µM) would overcome this issue. Our results showed that excess yDcpS completely removed the m7G cap from RNAs regardless of their linear or secondary structures (Fig. 2e). Above experiments established the careful design of the CompasSeq workflow, including the choice and order of decapping enzyme treatments. The workflow started with a cascade of enzymatic reactions, followed by adapter ligation (Fig. 1a). To test the efficacy of adapter ligation, we subjected 38 nt NAD-RNA and m7G-RNA to CompasSeq. As evidenced by the accumulation of an upper band (Fig. 2f), we demonstrated the specificity of cap-processing enzymes in combination with adapter ligation for NAD- and m7G-capped forms, respectively.The feasibility of CompasSeq analytical platformTo evaluate the feasibility of the CompasSeq analytical platform, encompassing both experimental and computational frameworks, we undertook several validation steps. First, we synthesized long RNA spike-ins (500 nt) with different sequences and cap structures, including NAD, FAD, dpCoA, m7G, and pppA. We mixed 20 µg of total RNA from mouse livers with 1 ng of each spike-in and subjected this mixture to the CompasSeq experiment, followed by qPCR analysis (Fig. 3a). Our results indicated that only NCIN, and not m7G- and ppp-RNAs, were enriched in the epitranscriptome-capturing group, whereas both NCIN- and m7G-RNAs were significantly enriched in the transcriptome-capturing group (Fig. 3b). Second, we aimed to measure the percentage of NCIN-capped RNA spike-ins. We excluded FAD-RNA spike-in from this and subsequent experiments due to its spontaneous decapping, which might confound calculations. Instead, we synthesized paired RNA spike-ins with identical sequences (500 nt) but different cap structures, including NCIN (NAD and dpCoA) or canonical m7G cap—each followed by a polyadenylated tail. It is likely that native transcripts exist in both NCIN- and m7G-capped versions, with proportions varying by gene. This experimental setup was intended to mimic endogenous transcripts that have 5% NCIN-capped modification (mixed with 95% m7G-RNA) or 0% NCIN-capped modification (i.e., 100% m7G-RNA) spike-in representing genes with no NCIN capping, thereby allowing us to assess the specificity of CompasSeq. We mixed 20 µg of total RNA from 2-month mouse livers with 0.014 ng of synthetic spike-in RNAs and subjected this mixture to the CompasSeq experiment, followed by polyA-selected RNA sequencing.Fig. 3: The feasibility of CompasSeq.a Schematic workflow of gene-specific assessment of NCIN modification with CompasSeq platform by qPCR. Total RNAs from 2-month-old mouse livers were mixed with polyadenylated spike-in RNAs (500 nt), followed by enzymatic reactions of CompasSeq for qPCR analysis. ppp-RNA (500 nt) was included as a baseline control. b qPCR analysis showed that synthetic RNAs (500 nt) containing NCIN caps, including NAD, FAD, and dpCoA, could be significantly and selectively enriched by the procedure designed for capturing NCIN-RNAs, whereas ppp- and m7G-RNA did not exhibit such enrichment (left). In the procedure designed for capturing m7G- and NCIN-RNAs, all RNAs, except ppp-RNA, were significantly enriched (right). Data were shown in mean ± s.e.m. (n = 3 independent experiments). Two-tailed Welch′s t-tests were used to infer p values. n.s. denoted not significant. c RNA-seq experiment of spike-in RNAs determined the specificity of CompasSeq. A total of three different spike-in RNAs with different capped forms were mixed with total RNAs from 2-month-old mouse liver, and the mixture of RNAs was subjected to CompasSeq. For each capped form, two polyadenylated spike-in RNAs with identical sequences (500 nt) but have either NCIN (NAD or dpCoA) or m7G cap, were used to mimic endogenous NCIN-capped transcript with different capped forms. As shown in the bar chart, NCIN-RNAs, but not m7G-RNA, were significantly and selectively enriched with NCIN capping levels close to their expected ratio. Data were shown in mean ± s.e.m. (n = 3 independent experiments). d RNA-seq experiment of spike-in RNAs determined the sensitivity of CompasSeq. A total of five different polyadenylated spike-in RNAs (500 nt) with different proportions of NAD-RNA were mixed with total RNAs from 2-month-old mouse livers, and the mixture of RNAs was subjected to CompasSeq. Two spike-in RNAs with identical sequences but have either NAD or m7G cap, were used to mimic endogenous NCIN-capped transcript with different NCIN capping levels. Scatter plot showing the concordance between expected and observed NCIN capping levels from spike-in RNAs. The Pearson correlation coefficient and the associated two-sided p value were shown in the top right corner. Data were shown in mean ± s.e.m. (n = 3 independent experiments). Exact p values and source data are provided as a Source Data file.Full size imageImportantly, we added pre-treated total RNA from Drosophila at a mass ratio of 50:1 for each sample, serving as a baseline control after adapter ligation, a critical procedure for downstream computational analysis (Supplementary Fig. 2a). We reasoned that exogenous Drosophila RNA spike-ins could be utilized to reflect and correct unwanted variations introduced by experimental procedures. In principle, Drosophila spike-ins, when added equally to each sample, should maintain a constant proportion across samples, capturing either NCIN- or all-capped forms. Given the fact that NCIN-capped RNAs constituted only a small proportion of the total RNA pool, read counts from the NCIN-capped library should be lower than those from the all-capped library. In practice, we observed that the proportion of reads from the Drosophila genome was higher in samples capturing NCIN-capped than all-capped forms (Supplementary Fig. 2b). Therefore, we employed a stepwise normalization procedure to address these variations, to enable the recovery of the theoretical distribution of NCIN-capped RNAs. We found that raw read counts displayed similar distributions between NCIN-capped RNA and all-capped RNA libraries (Supplementary Fig. 2c). Upon normalizations, gene measurements from all-capped forms exhibited higher values compared to those from NCIN-capped libraries, reflecting the intrinsic distribution of NCIN-capped transcripts (Supplementary Fig. 2c). Correlation analysis revealed that normalized gene measurements exhibited improved reproducibility compared to the raw read counts (Supplementary Fig. 2d). However, NCIN capping estimated from genes with relatively low expression exhibited higher variation compared to highly expressed genes, an observation potentially caused by genes with low read counts that tend to suffer from high variability40. To address this issue, we introduced the empirical Bayes shrinkage strategy. As shown in Supplementary Fig. 2e, CompassAnalyzer could significantly mitigate such variations. Comparatively, exogenous total RNA spike-ins added at the beginning step failed to correct PCR amplification biases, as shown by the poor performance of NCIN-capped RNA proportions (Supplementary Fig. 2f).Using CompassAnalyzer, we were able to analyze synthetic RNA spike-ins with known ratios, revealing that this platform not only captured metabolite-modified RNAs but also their relative ratios for the 5% of NCIN-capped spike-in (Fig. 3c). Background enrichment for the 100% m7G-RNA spike-in was significantly lower (less than 2%) compared to the 5% NCIN-capped spike-in (Fig. 3c), reflective of the specificity. We further employed spike-in RNAs with ascending ratios of NAD-capped forms: 0% (i.e., 100% m7G-capped), 1% (i.e., 99% m7G-capped), 5% (i.e., 95% m7G-capped), 10% (i.e., 90% m7G-capped), and 20% (i.e., 80% m7G-capped) (Supplementary Fig. 2g). The CompasSeq method, with its streamlined experimental and computational procedures, could reach coherent results between observed and expected NCIN capping levels (Pearson′s r = 0.998, Fig. 3d), thus illustrating the feasibility of this method in capturing and quantitatively measuring endogenous metabolite-modified RNAs.NCIN-RNA epitranscriptomeWith the CompasSeq analytical platform, we characterized NCIN-RNAs from young (2-month) and aged (18-month) mouse liver tissues. Following standard quality filtering procedures, each of the three independent biological replicates yielded an average of ~21.29 million high-quality reads that uniquely aligned to the reference genome (Supplementary Fig. 3a). PCA showed that the biological replicates of each sample type clustered together, indicating high reproducibility (Supplementary Fig. 3b). Dataset assessments corroborated that sequencing saturation has been reached (Supplementary Fig. 3c). Based on the background enrichment signal detected by 100% m7G-RNA spike-in, which was less than 2% (Supplementary Fig. 3d), we chose 2% ratio as the cutoff. This criterion led us to identify a total of 5680 highly reproducible NCIN-RNAs (Supplementary Data 2; Pearson′s r = 0.895 of the three biological replicates, Supplementary Fig. 3e), accounting for ~47.2% of the mRNA transcripts expressed in young mouse liver tissues. Besides RNA-seq, we also applied RT-qPCR. By analysis of three randomly selected genes, Ndufb10, Fkbp8, and Fads1, we calculated the relative abundance between NCIN-capped and all-capped groups (Supplementary Fig. 4a). The results obtained by RT-qPCR were in good congruence with those from high-throughput sequencing (Supplementary Fig. 4b), suggesting that CompasSeq can also be used for gene-specific assessment. Furthermore, comparison between the NCIN-RNAs, including, but not limited to, NAD-capped RNAs, identified in the current study with NAD-capped RNAs identified by our recently published ONE-seq method24 revealed a significant overlap (Supplementary Fig. 4c). For those NAD-RNAs uniquely identified by ONE-seq, their expression levels were significantly low (Supplementary Fig. 4d).We briefly analyzed the NCIN-capped epitranscriptome from young mouse livers. By computing NCIN capping levels, we observed that most genes producing NCIN caps had modification ratios between 5% and 30%, with some genes exhibiting exceptionally high levels, such as Rnr1 (91.1%) (Fig. 4a). NCIN capping predominantly occurred on protein-coding genes, but also extended to pseudogenes and noncoding RNAs (Fig. 4b). NCIN-RNAs were derived from genes located on autosomes, the X chromosome, and the mitochondrial genome, but not from the Y chromosome (Fig. 4c). Based on the extent of modification, we divided NCIN-RNAs into three categories: high (30–100%, e.g., Fads1), medium (12.5–30%, e.g., Fkbp8), and modest levels (2–12.5%, e.g., Ndufb10) (Fig. 4d and Supplementary Fig. 4e). To inspect biological function, we performed pathway enrichment analysis and found that genes with modest NCIN capping were mainly involved in mitochondria-related events and energy metabolism, whereas genes with medium and high NICN capping levels were enriched in macroautophagy, vesicle-related processes, and cellular responses (Fig. 4e and Supplementary Data 3).Fig. 4: The quantitative NCIN-capped epitranscriptome from mouse liver.a The distribution of NCIN capping levels of each gene from 2-month-old mouse livers (n = 5462). Box plot represented the median (center line), the first and third quartiles (box limits), and the whiskers extended to 1.5 times the interquartile range. Data points beyond the whiskers were shown as outliers. b Pie chart showing the gene types of identified NCIN-RNAs. c Bar chart showing the chromosome distribution of identified NCIN-RNAs. d The distribution of NCIN capping levels of genes with modest (2–12.5%), medium (12.5–30%), and high (30–100%) NCIN capping levels. e Pathway enrichment analysis of genes with modest, medium, and high NCIN capping levels. Numbers represented the count of NCIN-RNA enriched in each pathway. Multiple testing corrections were performed using gprofiler2’s built-in “g_SCS” method. The gray dashed line indicated the cutoff of 0.05 adjusted p value. Exact p values and source data are provided as a Source Data file.Full size imageValidation of CompasSeq by an orthogonal methodTo validate the CompasSeq platform by an orthogonal approach, we introduced two 5′−3′ exoribonucleases: dRai1 and Xrn1. Our rationale was based on the above observation that dRai1 specifically decaps NCIN-, but not m7G-capped RNAs, leading to 5′-to-3′ cleavage (Fig. 2b). This effect could be further enhanced by Xrn1 towards degradation, whereas m7G-capped transcripts would largely remain intact (Fig. 5a). We first validated the decapping activity of dRai1 on FAD- and dpCoA-RNAs (see Fig. 2b for NAD-RNA) by observing the appearance of lower bands on the UREA-PAGE gel electrophoresis (Fig. 5b). Second, we subjected 106 nt NAD-RNA and m7G-RNA spike-ins to dRai1 and Xrn1 treatment, which showed selective degradation of NAD-, but not m7G-capped, RNA (Fig. 5c). Prompted by this evidence, we hypothesized that the relative level of NCIN-capped transcripts could be determined by the reduction in abundance after treatment with dRai1 and Xrn1, compared to a mock treatment that includes both NCIN- and m7G-capped forms. We named this method the quantitative exoribonuclease reduction assay (Fig. 5d). We subjected 1 μg of total RNA from mouse livers mixed with 0.1 ng of synthetic RNA spike-ins (500 nt), including m7G, NAD, FAD, dpCoA, and pppA, to the assay. Consistently, only NCIN-RNAs exhibited significant degradation, while m7G-RNA remained unchanged (Fig. 5e). To assess whether the assay could measure NCIN capping levels, we mixed 1 μg of total RNA with 500 nt spike-in RNAs containing ascending ratios of NAD relative to m7G-capped forms: 0%, 10%, 30%, 50%, 70%, and 90%. We observed that relative capping levels measured by the assay corresponded to their expected levels (Fig. 5f). Furthermore, we successfully validated the NCIN capping levels of five randomly selected genes: Idh1, F13b, Ssr3, Rnr2, and Rnr1, as identified and quantified by CompasSeq (Fig. 5g, h). As a control experiment, we confirmed the absence of NCIN capping in five genes not identified by CompasSeq (Apom, Cox6a1, Dbi, Ddt, and Lrg1) (Supplementary Fig. 5). Taken together, our results using quantitative exoribonuclease reduction assay strongly supported the reliability of CompasSeq for the identification and ratio quantification of NCIN-RNAs.Fig. 5: Validation of CompasSeq by quantitative exoribonuclease reduction assay.a The validation strategy for the treatment with dRai1 and Xrn1. b UREA-PAGE gel showing products of processing of 38 nt FAD- and dpCoA-RNA 5′-ends by dRai1. c UREA-PAGE gel showing products of processing of 106 nt NAD- and m7G-RNA by dRai1 and Xrn1. d Schematic workflow of gene-specific assessment of NCIN modification with quantitative exoribonuclease reduction assay. e Synthetic RNAs (500 nt) containing NCIN caps (NAD, FAD, and dpCoA) could be significantly and selectively degraded by dRai1 and Xrn1, whereas m7G-RNA did not exhibit such degradation. Data were shown in mean ± s.e.m. (n = 3 independent experiments). Two-tailed Welch′s t-tests were used to infer p values. n.s. denoted not significant. f Synthetic RNAs (500 nt) with ascending ratios of NAD-capped forms: 0% (i.e., 100% m7G-capped), 10%, 30%, 50%, 70%, and 90% were used for the quantitative exoribonuclease reduction assay. The red and blue numbers represented the expected and observed NCIN capping levels, respectively. The observed NCIN capping levels were calculated as the difference between 1 and the relative levels from enzymatic treatment measured using the assay. Data were shown in mean ± s.e.m. (n = 3 independent experiments). g Five endogenous genes from mouse liver were assessed by the quantitative exoribonuclease reduction assay. Data were shown in mean ± s.e.m. (n = 3 samples). h Comparison between the result obtained by quantitative exoribonuclease reduction assay and CompasSeq. The red bars represented NCIN capping levels obtained from CompasSeq (left axis), while the blue bars illustrated NCIN capping levels calculated from the assay (right axis). Data were shown in mean ± s.e.m. (n = 3 samples). Two-tailed Welch′s t-tests were used to infer p values. n.s. denoted not significant. Exact p values and source data are provided as a Source Data file.Full size imageCharacterization of metabolite-capped RNA at the transcript resolutionWe proceeded to characterize newly identified NCIN-RNAs at the transcript resolution. Eukaryotic genomes feature transcript isoforms that can be derived from common gene loci. In-depth analysis revealed that, for most NCIN-RNAs (98.5%), their isoforms, if transcribed, all contained NCIN caps, e.g., Serpinf2, although the actual capping level varied among individual isoforms. Only 1.5% of NCIN-RNAs had select isoforms being modified, e.g., Zfp36 (Supplementary Fig. 6a, b and Supplementary Data 4). This observation prompted us to explore the link between NCIN capping and gene transcription, including transcript abundance, gene length, number of introns, alternative transcription start sites (ATSS)41, alternative splicing (AS)42, and alternative polyadenylation (APA)43. First, by dividing NCIN-RNAs into 10 deciles based on their NCIN capping levels, we observed that transcript abundance was negatively correlated with NCIN capping levels, while transcript length and number of introns were positively correlated (Supplementary Fig. 6c). Second, we found that RNAs with high levels of NCIN capping tended to possess ATSS (Supplementary Fig. 6d). Third, we noted that NCIN capping was negatively associated with AS events, including exon skipping, alternative 5′ donor and 3′ acceptor sites, mutually exclusive exons, and intron retention (Supplementary Fig. 6e). Fourth, we assessed the extent of APA by the percentage of distal poly(A) site usage index (PDUI). This analysis illustrated that transcripts containing NCIN caps were more likely to have APA compared to RNAs without NCIN modification. Consistently, NCIN capping levels positively correlated with distal poly(A) site usage (Supplementary Fig. 6f). NCIN capping represents an evolutionarily conserved RNA modification5,12,27,44,45. To explore the potential impact of sequence conservation, we utilized PhastCons scores, which measure conservation based on sequence similarity46, to calculate the conservation score of each isoform. Our findings indicated that transcripts with low sequence conservation exhibited less metabolite capping, whereas transcripts with high sequence conservation tended to be highly modified (Supplementary Fig. 7a). Prior research has indicated that proximal RNA secondary structures can affect RNA capping47. Therefore, we analyzed sequence features within 100 bp downstream of the transcription start site (TSS), assessing potential structural influences based on minimum free energy (MFE)48. Our analysis revealed that RNAs with NCIN caps tended to adopt predicted structures with lower MFE than those with m7G caps. Furthermore, RNAs with higher NCIN capping exhibited even lower MFE levels (Supplementary Fig. 7b, c).NCIN capping on transcript stability and translation efficiencyTo investigate how NCIN capping affects RNA stability, we generated Genome-wide Mapping of Uncapped and Cleaved Transcripts (GMUCT)49 datasets to assess mRNA degradation in mouse livers. All GMUCT libraries exhibited read enrichment at the 3′ end of mRNAs (Supplementary Fig. 7d), and GMUCT read abundance per transcript was positively correlated with mRNA levels quantified by RNA-seq (Supplementary Fig. 7e), which aligned with previously reported degradome sequencing patterns50. Using the GMUCT and RNA-seq data, we assessed transcript stability using the “proportion uncapped” metric (log2[GMUCT/RNA-seq]), where higher values suggest lower stability. We analyzed the stability of RNA transcripts with and without NAD24 and NCIN caps (CompasSeq data from the current study), respectively. Our analysis revealed that NAD-capped RNAs in mouse liver exhibited a higher proportion of 5′ degradation intermediates, suggesting greater instability compared to RNAs without NAD modification (Supplementary Fig. 7f). These observations were in line with findings in Arabidopsis, where NAD capping was found to be associated with reduced RNA stability51. In a sharp contrast, NCIN-capped RNAs from CompasSeq data showed a significantly lower proportion of degradation intermediates (Supplementary Fig. 7g). These data suggested that NCIN caps, as the compendium of NAD, FAD, dpCoA caps, etc., tended to occur at transcripts with higher stability. To examine the impact of NCIN capping on translational efficiency, we calculated ribosome occupancy (RO) using available Ribo-seq and RNA-seq data from mouse liver (GEO: GSE12398152). RO, an estimate of the translation efficiency, was computed as the log-ratio of Ribo-Seq and RNA-Seq read counts. Our analysis indicated that NCIN capping was positively correlated with translation efficiency (Supplementary Fig. 7h). Together, these findings suggested that NCIN caps were linked with transcripts exhibiting greater stability and higher translation efficiency.NCIN capping dynamics during agingThe number of genes harboring NCIN caps and their relative capping levels exhibited a significant decrease during aging (Fig. 6a). Moreover, we quantified the usage of metabolite caps by calculating the weighted sum of NCIN capping levels normalized to transcript abundance. Additionally, metabolite-capped RNA constituted ~8.3% of the mRNA transcriptome in young animals, and this proportion declined to 6.5% in aged mice, though the difference was not statistically significant (Fig. 6b).Fig. 6: Features of NCIN-capped epitranscriptome.a Violin plot showing the global dynamics of NCIN-RNAs during aging. Box plots represented the median (center line), the mean (red dots), the first and third quartiles (box limits), and the whiskers extended to 1.5 times the interquartile range. Data points beyond the whiskers were shown as outliers. Statistical analysis was performed using the nonparametric Mann–Whitney rank test (two-sided). b Pie chart showing the percentile usage of RNA caps (NCIN caps and non-NCIN cap) between young and old mouse livers. Statistical analysis was performed using a chi-squared test (p = 0.83, two-sided). c Scatterplot showing NCIN capping levels in young and old mouse livers. Differentially modified genes were defined as absolute fold change of NCIN capping levels ≥ 1.5 in old compared to young samples. d Scatterplot showing log2-transformed normalized gene expression levels in young and aged mouse livers. Differentially expressed genes were defined as absolute fold change of normalized read counts ≥ 1.5 and adjusted p values