Multiplexed direct RNA sequencing using RB‑dRNAseq with RNA barcodes

1 Nanopore direct RNA sequencing (dRNA-Seq) simultaneously resolves RNA 2 sequences, poly(A) tail lengths, and base modifications at the single-molecule level. 3 Here, we developed RB-dRNAseq, a multiplexed dRNA-Seq method using 4 RNA-barcoded reverse transcription adapters (RTAs) compatible with Oxford 5 Nanopore Technologies (ONT) workflows, enabling highly efficient sample 6 demultiplexing based solely on basecalled sequences. On the RNA004 platform, 7 RB-dRNAseq achieved >99.5% demultiplexing accuracy, with 81% of reads assigned 8 in 6-plex samples. This represents a 16%-24% relative increase in read recovery 9 compared to conventional DNA barcoding methods that rely on basecalling-based 10 demultiplexing. To evaluate the sensitivity of RB-dRNAseq, we performed 11 sequencing with low-input samples, demonstrating that even with 1 ng of total RNA, 12 gene expression (r=0.84) and m6A modification (r=0.91) profiles remained highly 13 consistent with high-input benchmarks. We further explored the feasibility of 14 RB-dRNAseq in ultra-low input scenarios by applying it to mouse zygotes and 2-cell 15 stage embryos. The method successfully captured dynamic shifts in gene isoform 16 expression and multiple RNA modifications (m6A, m5C, pseU, and Inosine) during 17 zygotic genome activation (ZGA), uncovering site-specific regulatory patterns. 18

19 The nanopore sequencing platform enables direct sequencing of native 20 RNA(Stark et al. 2019; Jain et al. 2022; Monzó et al. 2025), concurrently resolving 21 base sequence(Xin et al. 2021; Zhang et al. 2021), poly(A) tail length(Soneson et al. 22 2019; Parker et al. 2020; Roach et al. 2020), and RNA base modifications for 23 individual RNA molecules(Garalde et al. 2018; Begik et al. 2021; Lucas and Novoa 24 2023; Lucas et al. 2023; Liu and Conesa 2025). However, the nanopore dRNA-Seq 25

for low-input samples remain scarce, primarily due to several technical 26 challenges(Stark et al. 2019; Pardo-Palacios et al. 2024; Liu and Conesa 2025). First, 27 RNA extraction protocols for trace samples are underdeveloped. Second, current 28 library preparation and sequencing workflows demand high RNA 29 input(Pardo-Palacios et al. 2024). For instance, the recommended input for 30 dRNA-Seq using the ONT RNA004 flow cell is 1 μ g total RNA or 300 ng poly(A)+ 31 RNA, far exceeding the RNA content of most low-input samples and limiting its 32 application. 33 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint Multiplexing samples using barcoded Reverse Transcription Adapters (RTAs) 34 offers a viable strategy to increase the effective starting material and thus ensure 35 sufficient sequencing yield in dRNA/i2 Seq(Foord et al. 2023). However, ONT has not 36 yet officially released a barcoding method for dRNA-Seq. A key challenge is that the 37 RTA sequence is DNA, and DNA barcodes are co-sequenced with the RNA strand, 38 making them difficult for RNA basecalling software to recognize, thus hindering their 39 use in dRNA-Seq. Recent academic developments include several multiplexed 40 dRNA-Seq techniques using DNA-barcoded RTAs. For example, Smith et al. (2020) 41 proposed DeePlexiCon, which designs a 20 bp DNA barcode within the RTA and 42 utilizes deep learning models trained on the current signal of the barcode sequences 43 for demultiplexing. On the RNA002 platform, DeePlexiCon classified 60% of reads 44 with 99.9% accuracy using 4 barcodes(Smith et al. 2020). Subsequently, Toorn et al. 45 released WarpDemuX, employing a similar DNA-barcoded RTA scheme and 46 extending the current signal-based deep learning demultiplexing model to 12 barcodes, 47 further increasing the sample throughput(van der Toorn et al. 2025). Additionally, 48 SeqTagger improved upon DeePlexiCon by expanding barcode diversity and adapting 49 the approach to the RNA004 platform, and it performs demultiplexing via basecalling 50 with a custom-trained Bonito model(Pryszcz et al. 2025). While DeePlexiCon, 51 WarpDemuX, and SeqTagger provide feasible barcoded dRNA-Seq solutions and 52 allow researchers to use pre-trained models, their flexibility is limited. Customization 53 requires additional deep learning model training, increasing the barrier to use. 54 Furthermore, these methods necessitate retraining following major ONT upgrades to 55 reagents, flow cells, or basecalling algorithms. 56 Here, we developed RB-dRNAseq, which employs RNA-barcoded RTAs 57 compatible with the ONT dRNA-Seq library prep, enabling high-precision 58 demultiplexing post-sequencing based solely on the basecalled sequence. We 59 demonstrated that the RNA barcoding strategy significantly increased the 60 demultiplexing efficiency compared to DNA barcoding, and the performance 61 improved from RNA002 to RNA004 platform. We show that RB-dRNAseq accurately 62 capture the isoform expression, poly(A) tail length and epitranscriptome information 63 across different total RNA input amounts. Finally, we applied RB-dRNAseq to mouse 64 zygote and 2-cell stage embryos, resolving isoform expression and modification 65 dynamics at the single-molecule level during zygotic genome activation (ZGA). 66

67 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint High accuracy and specificity of demultiplexing DRS data using the RNA 68 barcoding strategy 69 For RNA-barcoded RTA design in RB-dRNAseq, we developed OligoA 70 (5'-10DNA-24RNA-3DNA-3') where the 20 RNA barcode bases are flanked by DNA 71 sequences, and OligoB is entirely DNA (Fig. 1A). OligoA and OligoB are synthesized 72 separately and annealed to form the complete barcoded RTA (Supplementary Table 1). 73 This single-stranded RNA barcode design reduces synthesis costs and the three DNA 74 bases at the 3' end of OligoA ensure compatibility with the ONT ligation step and 75 mitigate RNA degradation risk. We replaced the standard RTA component with our 76 RB-barcoded RTA, all other library preparation and sequencing steps strictly followed 77 the ONT protocol when using total RNA input (Fig. 1B). 78 We first tested the RNA barcode strategy on RNA002 and this largely improved 79 the reads demultiplexing ratio (70.44%) compared with the DNA barcodes (54.51%, 80 Fig. 1C). The ratio of demultiplexed reads further improved on RNA004 (87.77%) 81 using the RNA barcode strategy (Fig. 1C, Supplementary Table 2). We examined the 82 impact of barcode alignment e-value thresholds on demultiplexing rate and accuracy. 83 As the e-value threshold became stricter, the accuracy changed minimally 84 (99.14%-99.88%), but the demultiplexing rate gradually decreased (Fig. 1D, 85 Supplementary Fig. 1A-B, Supplementary Table 2). We finally selected an e-value 86 threshold of 0.05 for downstream analysis. 87 We further evaluated the impact of different ONT Dorado basecalling modes on 88 demultiplexing rate and accuracy. The Dorado sup mode significantly increased the 89 barcode demultiplexing rate compared to the fast mode, which the accuracy slightly 90 decreased (Fig. 1E, Supplementary Fig. 1C-D, Supplementary Table 2). Additionally, 91 increasing the number of barcode types from 2 to 6 also caused the demultiplexing 92 rate decreased from 89.2% to 81.4% (Fig. 1F, Supplementary Fig. 1E-F, 93 Supplementary Table 2). We checked whether the failed reads (without barcode 94 sequence recognition) were caused by the sequencing quality or the read length, but 95 neither of the factors exhibited differences between barcoded reads and non-barcoded 96 reads (Fig. 1G-H). This might be caused by defects of the current basecalling 97 algorithm at the barcode region. Future optimization of the basecalling algorithm may 98 improve the overall demultiplexing rate, supporting the use of more barcode types. 99 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint 100 Figure 1. Schematic and performance evaluation of RB-dRNAseq. (A) Schematic 101 of barcoded RTA design. (B) RB-dRNAseq library preparation workflow. (C) Bar plot 102 comparing demultiplexing rates of RNA002 DNA barcode, RNA002 RNA barcode, 103 and RNA004 RNA barcode strategies (hac basecalling mode, e-value<0.05). (D) Line 104 plot (red, demultiplexing rate) and bar plot (blue, accuracy) for RNA004 RNA 105 barcode under different e-value thresholds (super basecalling mode). (E) Line plots 106 comparing demultiplexing rate and accuracy for RNA004 RNA barcode under 107 different Dorado basecalling modes (e-value<0.05). (F) Line plots comparing 108 demultiplexing rate and accuracy for RNA004 RNA barcode with different numbers 109 of barcodes (sup basecalling mode, e-value<0.05). (G) Read length distributions for 110 reads with demultiplexed barcodes vs. undemultiplexed reads (sup basecalling mode, 111 e-value<0.05). (H) Sequencing quality (Q-score) distributions for reads with 112 demultiplexed barcodes vs. undemultip lexed reads (sup basecalling mode, 113 e-value<0.05). 114 Validation of RB-dRNAseq across different total RNA input 115 To evaluate the performance of RB-dRNAseq across varying RNA input 116 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint amounts from different species, we extracted total RNA from human HEK293T cells 117 and mouse embryonic stem cells (mESCs), respectively. We used 1 ng, 10 ng, 100 ng, 118 500 ng of total RNA as the starting material, corresponding to four different barcodes 119 (Fig. 2A). Barcode demultiplexing demonstrated that RB-dRNAseq accurately 120 recapitulated the output read counts corresponding to different input amounts (Fig. 2B, 121 Supplementary Table 3). Quality control analysis showed stable read length 122 distributions across input levels (Fig. 2C), and there were no significant difference 123 between the average poly(A) tail lengths detected in the low-input and high-input 124 RNA samples (Fig. 2D). 125 Then we compared the performance in gene and transcript detection by 126 RB-dRNAseq when starting with different amounts of input RNA. Both the number 127 of genes and isoforms largely increased from 1ng to 500ng total RNA as input 128 (Supplementary Fig. 2A-B). This should attribute to lower sample input leading to 129 fewer sequencing output reads. We confirmed that the number of detected genes 130 largely increased with the number of obtained reads under current shallow sequencing 131 depth (Supplementary Fig. 2C-D). For 1ng total RNA, we detected 1828 genes, 1371 132 transcripts in human and 2531 genes, 1791 transcripts in mouse (Supplementary Table 133 3). Although the number of genes and isoforms detected by 1ng total RNA is only 134 one-tenth of that detected by 500ng total RNA, the correlation between the gene 135 expressions obtained from the two groups is still high (correlation coefficient value 136 0.81-0.86, Fig. 2E-F, Supplementary Fig. 3A-B). This indicates reliable transcript 137 quantification in low-input samples using RB-dRNAseq. 138 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint 139 Figure 2. Performance of RB-dRNAseq with different total RNA input. (A) 140 Experimental design for multiplexed RB-dRNAseq with varying input amounts. (B) 141 Line plot comparing sequencing read counts across input amounts. (C) Boxplot of 142 read length distributions across input amounts. (D) Boxplot of poly(A) tail length 143 distributions across input amounts. (E-F) Scatter plot showing gene expression (E) 144 and isoform expression (F) correlation between different HEK293T total RNA input 145 indicated by the RTA barcodes. 146 As RB-dRNAseq also reports RNA modifications at single-base resolution, we 147 evaluated the accuracy of m6A detection in each group. Similar to expression 148 detection, the number of detected m6A sites highly relies on the amount of input 149 sample (Supplementary Fig. 4A-B). Metagene distribution profiles of m6A sites all 150 showed enrichment in the 3’ UTR, consistent to the classical m6A distribution (Fig. 151 3A-B). However, the 100 and 500ng groups displayed higher density of m6A sites in 152 the 3’UTR than the 10 and 1ng groups, indicating that the low-input samples would 153 prefer lose the capture of m6A in this region. Interestingly, the quantification of the 154 m6A modification level is not significantly affected by the initial amount of the input 155 RNA, as the correlation values among different groups of samples did not changed 156 obviously (Fig. 3C-D). These data robustly demonstrates accurate capture of m6A 157 sites and precise quantification of the modification levels in low-input samples using 158 RB-dRNAseq. 159 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint 160 Figure 3. Comparison of m6A detection by RB-dRNAseq across different total 161 RNA input. (A-B) Metagene distribution of m6A sites for HEK293T (A) and mESC 162 (B) across different total RNA input indicated by the RTA barcodes. (C-D) Scatter 163 plot showing the correlation of modification level of each m6A site between different 164 HEK293T (C) and mESC (D) total RNA input indicated by the RTA barcodes. 165 RB-dRNAseq in low-input cell samples 166 dRNA-Seq has been only applied for large amount of total RNA or poly(A) RNA 167 in previous studies. We proved that RB-dRNAseq accurately capture expression and 168 modification of as few as 1ng total RNA. However, in practical applications, such 169 amount of samples are generally not subjected to RNA extraction. So added cell lysis 170 and RNA capture step before ligating RTAs and applied the workflow in mouse 171 zygote and 2-cell embryos. We started with as few as 1 embryos in each sample 172 replicate. Meanwhile, we added spike-in RNA in the library to maintain the activity of 173 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint the sequencing holes and ensure sufficient sequencing yield (Fig. 4A, Supplementary 174 Table 4). 175 Due to too few reads generated from several batches of samples, we merged 176 reads from the same stages and compared the transcript lengths and poly(A) tail 177 lengths at the zygote and 2-cell stages. The 2-cell embryos contains longer RNAs with 178 longer poly(A) tails (Fig. 4B-C), which is consistent to previous studies (e.g., Liu et 179 al.(Liu et al. 2023)). Further mapping revealed 2218 genes and 1405 RNA isoforms in 180 zygotes, and 3401 genes and 1759 isoforms in 2-cell embryos (Fig. 4D). Integration 181 with our previously published HIT-scISOseq data showed higher concordance in gene 182 expression between RB-dRNAseq and HIT-scISOseq for the same stages of samples 183 (Fig. 4E-G), validating the reliability of RB-dRNAseq data for expression 184 quantification of low-input cell sample. 185 186 Figure 4. RB-dRNAseq results for mouse zygote and 2-cell samples. (A) 187 Schematic of RB-dRNAseq for mouse zygote and 2-cell embryos. (B-C) Distribution 188 of read lengths (B) and poly(A) tail lengths (C) for zygote and 2-cell embryos. (D) 189 Bar plot showing number of genes and isoforms detected by RB-dRNAseq in zygote 190 and 2-cell embryos. (E-F) Scatter plot of gene expression correlation between zygote 191 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint and 2-cell samples using RB-dRNAseq (E) and HIT-scISO-seq (F). (G) Violin plot 192 showing distribution of gene expression correlation values between RB-dRNAseq and 193 HIT-scISOseq at the single-cell level for zygote and 2-cell samples. 194 Isoform expression and RNA modification changes during ZGA 195 The fact that ZGA serves as the foundation for the emergence of life has received 196 widespread attention. Previous single-cell methods enable us to investigate gene 197 expression changes during this process. Due to limitations in the current single-cell 198 epitranscriptome detection methods, little is known about how RNA modification 199 changes during ZGA. With RB-dRNAseq from zygote and 2-cell embryos, we 200 focused on differences in multiple RNA modifications at RNA isoform resolution. We 201 performed joint detection of multiple RNA modifications (m6A, m5C, pseU, Inosine) 202 on single molecules using ONT's Dorado basecaller (see Methods). Subsequently, the 203 ONT modkit tool was used to map read-level methylation information to RNA 204 isoforms, yielding single-base resolution modification information for each gene 205 isoform. For all four kinds of modifications analyzed, the modified sites in zygotes 206 showed higher modification levels than those in 2-cell embryos (Fig. 5A), consistent 207 with previous reports (e.g., Wu et al.(Wu et al. 2022)). Further comparison of 208 differentially modified sites (DMS, difference exceeds 5%) returned more sites 209 exhibiting reduced modification levels in 2-cell embryos. Specifically, a total of 2796 210 DMS were identified, and the number of DMS for m6A and pseU was significantly 211 higher than for m5C and Inosine (Fig. 5B). Interestingly, we observed different 212 modification changes for different sites on the same transcript. Take the 213 ENSMUST00000114431 from Btg4 as example, m6A and Inosine levels were 214 generally low on this transcript in both group of samples; the zygote showed higher 215 pseU level at pos167 than the 2-cell embryo, while the opposite trend was observed at 216 pos437 and pos520; the 2-cell embryo exhibited increased m5C level at pos517 and 217 decreased at pos648 (Fig. 5C). Importantly, we can determine how same and different 218 modifications are coordinated on the same RNA molecule (Fig. 5D). These results 219 fully demonstrate RB-dRNAseq's capability to resolve heterogeneity of multiple RNA 220 modifications across samples at isoform and single-base resolution. 221 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint 222 Figure 5. RNA modification analysis at single molecular and single base 223 resolution. (A) Density distribution of modification levels of detected sites in the 224 isoforms. (B) V olcano plots showing DMS (2-cell vs zygote) for different 225 modification types. (C) Methylation abundance profiles for different modification 226 types across the Btg4 transcript ( ENSMUST00000114431) in zygote and 2-cell 227 samples. (D) Single-molecule read visualization of co-occurrence of different 228 modification types at specific sites on the Btg4 transcript ( ENSMUST00000114431) 229 for zygote and 2-cell samples. 230 231 232 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint

233 Multiplex sampling has been the major challenge for dRNA-Seq. We developed 234 RB-dRNAseq, which addresses the core limitation of DNA barcodes in dRNA-Seq - 235 unrecognizability by basecalling software, by designing the barcode sequence as RNA. 236 Meanwhile, we designed the DNA-RNA chimera sequence to reduce costs while 237 ensuring compatibility with ONT workflows via the 3' DNA anchor. Compared to 238 DNA barcoding methods relying on current signal deep learning (e.g., DeePlexiCon), 239 RB-dRNAseq enables demultiplexing through simple sequence alignment, 240 eliminating the complexity of model training and platform adaptation, while 241 significantly improving accuracy (>99.5%) and multiplexing capacity (>81% 242 demultiplexing rate at 6-plex). RB-dRNAseq maintains stable transcript quantification 243 (gene expression correlation r=0.84) and modification detection (m6A abundance 244 correlation r=0.91) even at the ultra-low input of 1 ng, with no significant bias in read 245 length or poly(A) tail distribution. Combined with a direct post-lysis library 246 construction strategy for single cells/embryos, the method was successfully applied to 247 study mouse early embryonic development. The detected transcript expression 248 patterns during ZGA were highly consistent with HIT-scISOseq data, validating data 249 reliability. In this way, we are able to integrate multidimensional information, 250 including isoform expression, poly(A) tail length, and multiple modifications 251 concurrently on single molecules, revealing coordinated changes during the ZGA 252 process (e.g., on Btg4). 253 However, there are also some limitations for the current version of RB-dRNAseq. 254 First, the 24 nt single-stranded RNA barcodes used in RB-dRNAseq have relatively 255 high synthesis costs. Given the current 6-plex RNA barcode edit distances of 12-18 256 (Supplementary Table 5), future optimization could reduce barcode length, thereby 257 lowering costs. Additionally, the demultiplexing ratio could be improved, especially 258 for more sample numbers. The DNA sequences flanking the RNA barcode may cause 259 basecalling errors. Future optimization of algorithms or targeted model training is 260 needed to improve demultiplexing rates. Second, the current sequencing yields for 261 low-input samples (including both total RNA and cells) are quite low. We tried to add 262 high amounts of spike-in RNA to ensure the activity of the sequencing holes, thus to 263 generate more total reads. But the proportion of reads from target samples is relatively 264 low. The best way to solve this problem is to enhance the efficiency of ONT 265 sequencing. This may extend the application of RB-dRNAseq in single-cell and 266 spatial transcriptome analysis, and rapid pathogen detection, etc. 267 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint

268 Animals and single blastomere collection 269 We used 8 to 12 week-old C57BL/6J female mice and DBA/2NCrl male mice in 270 the experiment. The female mice were first injected with 7.5 IU of pregnant mare’s 271 serum gonadotropin (PMSG) (Ningbo SanSheng Biological Technology, Cat. 272 110044564), and with 7.5 IU of human chorionic gonadotropin (hCG) (Ningbo 273 SanSheng Biological Technology, Cat. 50030248) after 48 hours. After mating, the 274 embryos of each stage were collected at defined time periods after hCG injection: 22 275 to 24 hours (zygote), 46 to 48 hours (late 2cell). All animal experiments were 276 performed according to the guidelines of the Guangzhou national laboratory 277 (Guangzhou, China). Collection of single blastomeres at each stage was carried out as 278 previously described(Koo et al. 2020). 279 mRNA capture for total RNA and cell lysate 280 Total RNA was extracted from HEK293T and mESC (derived from C57BL/6J 281 mice) by RNA extraction kit (Yeasen, Cat. 19221ES ） respectively. For early mouse 282 embryos, including zygotes and 2 cells, single embryo was taken and placed into the 283 cell lysis buffer (0.8% Triton X-100 (Sigma-Aldrich, T9284), 2 mM DTT 284 (Sigma-Aldrich, 43815), 2U/µl of RNase inhibitor (TaKaRa, 2313A)). V ortex 285 vigorously for 10 seconds to ensure complete lysis. 286 Resuspend the mRNA capture magnetic beads (Vazyme, Cat. N401), take 10 μ L 287 of magnetic beads into a 200 μ L volume centrifuge tube, place it on a magnetic stand 288 for 2 minutes, and discard the supernatant. Add 40 μ L of 1× HB buffer (3× 289 Saline-Sodium Citrate Buffer (Sigma-Aldrich, S6639), 0.05% Tween-20 290 (Sigma-Aldrich, 655204), 2U/µl of RNase inhibitor) to wash the magnetic beads and 291 repeat the washing twice. Add 40 μ L of 2× HB buffer (6×SSC, 0.05% Tween-20, 292 2U/µl of RNase inhibitor) to resuspend the magnetic beads. 293 Add an equal volume of the magnetic beads in 2× HB buffer to total RNA or cell 294 lysate, mix gently, and then instantly centrifuge at low speed and place it in PCR 295 instrument, incubate at 65°C for 5 minutes, then place it on a rotating rack at room 296 temperature for 15 minutes. Leave it on a magnetic stand for 2 minutes and discard 297 the supernatant. Add 40 μ L of pre-cooled 1× HB buffer to resuspend the magnetic 298 beads, centrifuge at low speed instantaneously, and then place it on a magnetic stand 299 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint for 2 minutes, discard the supernatant. Repeat washing twice. Add 40 μ L of 300 pre-cooled 1× phosphate-buffered saline (PBS, Invitrogen, AM9625) to wash the 301 magnetic beads once. Remove the supernatant thoroughly. Add 10 μ L of RNase-free 302 water, flick and resuspend thoroughly, instantly centrifuge at low speed, place it in a 303 PCR instrument, incubate for 2 minutes at 80°C, quickly incubate on ice for 2 minutes, 304 then place it on a magnetic stand for 5 minutes, transfer the supernatant to a new 200 305 μ L volume centrifuge tube. 306 Design and assembly of custom reverse transcription adapters (RTAs) containing 307 Barcode 308 RTA is made of the forward chain RTA-F and the reverse chain RTA-R. The 309 forward strand RTA-F is a DNA-RNA chimeric strand composed of DNA bases and 310 RNA bases. Its structure is as follows: 5’-TCCN(24)TAGTAGGTTC-3’, where N(24) 311 is a 24nt-length RNA sequence, which is the barcode region of RTA. The other bases 312 are DNA bases. RTA-R consists of DNA bases and has a structure as follows: 313 5’-GAGGCGAGCGGTCAATTTTM(24)GGATTTTTTTTTTT-3’, where M(24)GGA 314 is the reverse complementary sequence of TCCN(24) in RTA-F. The Barcode 315 sequences are shown in supplemental table1. The RTA-F and RTA-R oligo are 316 synthesized by Sangon Biotech (Shanghai). 317 The RTA-F and RTA-R of each RTA were diluted to 10 μ M and 9 μ M with 318 annealing buffer (10 mM Tris-HCl, pH 7.5, 50 mM NaCl, 0.5 mM EDTA), 319 respectively. Take 5 μ L each and mix it and then denature and anneale. The reason for 320 using excessive forward strands in RTA assembly is that excessive forward strands 321 can ensure complete consumption of reverse strands, avoiding the reduction of 322 RNA-RTA products after annealing, resulting in competition with RTA during 323 subsequent reverse transcription. The assembly of custom RTA was carried out 324 according to the following PCR procedure: denaturation at 98°C for 3 minutes, 325 decreased by 1°C per minute until 25/i2 , and maintains at 4/i2 . RTA with different RNA 326 barcodes could be used for direct RNA sequencing library construction. 327 Sample multiplexing and library construction 328 According to the ONT Direct RNA Sequencing Product Manual (SQK RNA004), 329 the captured mRNA of each sample is ligated with RTA containing certain barcode, 330 and undergoes reverse transcription reaction to form an RNA-DNA hybrid strand. The 331 RNA/DNA hybrid is then ligated with RNA ligation adapter (RLA) to form a 332 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint complete sequencing structure. 333 Brief steps are as follows: For each sample, an RTA containing a specific 334 barcode is ligated with its RNA using NEBNext rapid ligase (NEB M0202M) at room 335 temperature for 30 minutes. Reverse transcription was then performed using 336 SuperScriptIV reverse transcriptase (ThermoFisher, 18090050). The prouducts of 337 reverse transcription with different RTAs were mixed together and then purified using 338 2× RNA Clean XP beads (Beckman, A63987). The purified product was ligated with 339 the RLA linker using NEB Nextflash ligase (NEB M0202M) at room temperature for 340 30 minutes. The reaction mixture was purified using 0.6× RNA Clean XP beads 341 (Beckman, A63987), washed twice with buffer (WSB), and the product was eluted 342 with elution buffer (EB). Then prepare the priming mixture according to the 343 instructions (SQK RNA004) and adjust it according to the volume of the library 344 adapted to the PromethION sequencer. Load the mixed solution on the machine into 345 the ONT sequencing cell that passes the quality check, set the sequencing parameters 346 in the MinKNOW software, and the sequencing lasts for more than 12 hours to obtain 347 sufficient sequencing data. 348 Basecalling and demultiplexing 349 After obtaining barcode-labeled DRS data, use ONT's official Dorado software 350 and set the "--no-trim" parameter to retain the sequencing adapter, and perform 351 basecalling on the original pod5 data. Three precision modes could be chosen: fast 352 (fast), high precision (hac), and super-high precision (sup). Next, the blastn program 353 was used to align the RTA barcode sequences to the DRS reads after basecalling. The 354 alignment results are filtered based on the alignment length greater than 9, the 355 alignment position (s.end) at the 5' terminal of barcode greater than 5, and the e-value 356 of the alignment. Finally, for each DRS read, the barcode with the smallest e-value in 357 alignment was selected as the barcode of the read. 358 Data availability 359 All sequencing data will be publicly released upon acceptance of the manuscript 360 for publication. Source data for Supplementary Table 2 are provided as a 361 supplementary data file. 362

363 X.-Y .F. was supported by grants from the Guangdong Provincial Pearl River 364 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint Talents Program (2021QN02Y747), Guangzhou science and technology elite "pilot" 365 project (SL2024A04J01788) and the Major Project of Guangzhou National 366 Laboratory (GZNL2023A02003, GZNL2024A03001). Z.-X.S. was supported by 367 grants from the National Natural Science Foundation of China (42107148). 368 Author Contributions Statement 369 X.-Y .F., and F.Z. conceived and designed the project; Z.-X.S., X.-Y .F., and 370 Q.-P.H. developed the experimental technology; Q.-P.H performed sequencing 371 experiments. Z.-X.S., and Q.-P.H. perform ed the informatics analysis; Z.-X.S., 372 Q.-P.H., and Y .-F.Z. coordinated data release and assisted with executing the pipeline. 373 Z.-X.S., Q.-P.H., and X.-Y .F. wrote the manuscript and created the figures; Y .-F.Z., 374 and F.Z. advised the study and revised the manuscript. All authors have read and 375 approved the final version of this manuscript. 376 Competing Interests Statement 377 X.-Y .F., Z.-X.S., and Q.-P.H. have filed a patent on the RB-dRNAseq method. 378

379 Begik O, Lucas MC, Pryszcz LP, Ramirez JM, Medina R, Milenkovic I, Cruciani S, 380 Liu H, Vieira HGS, Sas-Chen A et al. 2021. Quantitative profiling of 381 pseudouridylation dynamics in native RNAs with nanopore sequencing. 382 Nature Biotechnology 39: 1278-1291. 383 Foord C, Hsu J, Jarroux J, Hu W, Belchikov N, Pollard S, He Y , Joglekar A, Tilgner 384 HU. 2023. The variables on RNA molecules: concert or cacophony? Answers 385 in long-read sequencing. Nature Methods 20: 20-24. 386 Garalde DR, Snell EA, Jachimowicz D, Sipos B, Lloyd JH, Bruce M, Pantic N, 387 Admassu T, James P, Warland A et al. 2018. Highly parallel direct RNA 388 sequencing on an array of nanopores. Nature Methods 15: 201-206. 389 Jain M, Abu-Shumays R, Olsen HE, Akeson M. 2022. Advances in nanopore direct 390 RNA sequencing. Nature Methods 19: 1160-1164. 391 Koo B-K, Fan X, Tang D, Liao Y , Li P, Zhang Y , Wang M, Liang F, Wang X, Gao Y et 392 al. 2020. Single-cell RNA-seq analysis of mouse preimplantation embryos by 393 third-generation sequencing. PLOS Biology 18. 394 Liu T, Conesa A. 2025. Profiling the epigenome using long-read sequencing. Nature 395 Genetics 57: 27-41. 396 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint Liu Y , Zhao H, Shao F, Zhang Y , Nie H, Zhang J, Li C, Hou Z, Chen Z-J, Wang J et al. 397 2023. Remodeling of maternal mRNA through poly(A) tail orchestrates human 398 oocyte-to-embryo transition. Nature Structural & Molecular Biology 30: 399 200-215. 400 Lucas MC, Novoa EM. 2023. Long-read sequencing in the era of epigenomics and 401 epitranscriptomics. Nature Methods 20: 25-29. 402 Lucas MC, Pryszcz LP, Medina R, Milenkovic I, Camacho N, Marchand V , Motorin Y , 403 Ribas de Pouplana L, Novoa EM. 2023. Quantitative analysis of tRNA 404 abundance and modifications by nanopore RNA sequencing. Nature 405 Biotechnology 42: 72-86. 406 Monzó C, Liu T, Conesa A. 2025. Transcriptomics in the era of long-read sequencing. 407 Nature Reviews Genetics doi:10.1038/s41576-025-00828-z. 408 Pardo-Palacios FJ, Wang D, Reese F, Diekhans M, Carbonell-Sala S, Williams B, 409 Loveland JE, De María M, Adams MS, Balderrama-Gutierrez G et al. 2024. 410 Systematic assessment of long-read RNA-seq methods for transcript 411 identification and quantification. Nature Methods 21: 1349-1363. 412 Parker MT, Knop K, Sherwood A V , Schurch NJ, Mackinnon K, Gould PD, Hall AJW, 413 Barton GJ, Simpson GG. 2020. Nanopore direct RNA sequencing maps the 414 complexity of Arabidopsis mRNA processing and m6A modification. eLife 9. 415 Pryszcz LP, Diensthuber G, Llovera L, Medina R, Delgado-Tejedor A, Cozzuto L, 416 Ponomarenko J, Novoa EM. 2025. Rapid and accurate demultiplexing of 417 direct RNA nanopore sequencing data with SeqTagger. Genome Research 35: 418 956-966. 419 Roach NP, Sadowski N, Alessi AF, Timp W, Taylor J, Kim JK. 2020. The full-length 420 transcriptome of C. elegans using direct RNA sequencing. Genome Research 421 30: 299-312. 422 Smith MA, Ersavas T, Ferguson JM, Liu H, Lucas MC, Begik O, Bojarski L, Barton 423 K, Novoa EM. 2020. Molecular barcoding of native RNAs using nanopore 424 sequencing and deep learning. Genome Research 30: 1345-1353. 425 Soneson C, Yao Y , Bratus-Neuenschwander A, Patrignani A, Robinson MD, Hussain 426 S. 2019. A comprehensive examination of Nanopore native RNA sequencing 427 for characterization of complex transcriptomes. Nature Communications 10. 428 Stark R, Grzelak M, Hadfield J. 2019. RNA sequencing: the teenage years. Nature 429 Reviews Genetics 20: 631-656. 430 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint van der Toorn W, Bohn P , Liu-Wei W, Olguin-Nava M, Gribling-Burrer A-S, Smyth 431 RP, von Kleist M. 2025. Demultiplexing and barcode-specific adaptive 432 sampling for nanopore direct RNA sequencing. Nature Communications 16. 433 W u Y , Xu X, Qi M, Chen C, Li M, Y an R, Kou X, Zhao Y , Liu W , Li Y et al. 2022. 434 N6-methyladenosine regulates maternal RNA maintenance in oocytes and 435 timely RNA decay during mouse maternal-to-zygotic transition. Nature Cell 436 Biology 24: 917-927. 437 Xin R, Gao Y , Gao Y , Wang R, Kadash-Edmondson KE, Liu B, Wang Y , Lin L, Xing 438 Y . 2021. isoCirc catalogs full-length circular RNA isoforms in human 439 transcriptomes. Nature Communications 12. 440 Zhang J, Hou L, Zuo Z, Ji P, Zhang X, Xue Y , Zhao F. 2021. Comprehensive profiling 441 of circular RNAs with nanopore sequencing and CIRI-long. Nature 442 Biotechnology 39: 836-845. 443 (which was not certified by peer review) is the author/funder. All rights reserved. No reuse allowed without permission. The copyright holder for this preprintthis version posted May 17, 2026. ; https://doi.org/10.64898/2026.05.14.725048doi: bioRxiv preprint