To mitigate the challenge of rapid polynucleotide translocation through nanopores that complicates signal discrimination, we engineered a Pif1-like helicase through random amino acid mutations. Specifically, cysteine residues and unnatural amino acids were introduced into critical structural domains, including the tower (residues E264-P278, N296-A394), pin (K89-E105), and 1A domains (M1-L88, M106-V181) (Figure 1A). These modifications enhanced helicase-DNA binding and stabilized translocation at 400 bp/s, a rate comparable to ONT R9.4 flow cells (450 bp/s) and superior to earlier R9 iterations (250 bp/s) (Wang et al., 2021). The engineered helicase thus effectively regulates nucleic acid movement, improving the signal-to-noise ratio in downstream measurements.

In parallel, we redesigned the nanopore protein, derived from the CsgG/HfaB family of Pseudomonas species, by mutating amino acids 69–76 (KPTPASSF) to RPSPASAQ (Figure 1B). This variant exhibited improved structural rigidity and stronger nucleic acid affinity, leading to enhanced electrical signal stability. Together, the modified helicase and pore enabled sequencing with base calling accuracies approaching 90% and methylation detection correlations ≥96% with bisulfite sequencing (BS-seq) (see in following sections), comparable to ONT’s R9.4 platform (Wang et al., 2021).

To evaluate signal characteristics, we analyzed representative raw electrical traces from sequencing reads (Figure 1C). Signals outside the DNA-bound state exhibited higher current levels, while DNA-associated signals fluctuated within the 100–200 pA range. Zoomed-in views revealed sporadic drops, indicative of local signal instability—suggesting room for further pore optimization.

Next, we computed the average signal values across 100 intervals of 15-base segments for each canonical base (A, C, G, and T), and visualized the resulting distributions (Figure 1D). To explore the underlying structure of these signal patterns, we applied t-SNE dimensionality reduction (van der Maaten and Hinton, 2008) to over 400 normalized current signal data points. As shown in Figure 1E, signals associated with different bases formed partially separable clusters, indicating that some base-specific features were preserved. However, due to substantial signal noise and overlap among distributions, traditional statistical methods proved inadequate for robust base identification.

We extended this analysis to methylation-related signals by comparing CpG fully methylated (YF6418; all cytosines in CpGs as 5 mC) and fully unmethylated (YF6419) datasets from D. melanogaster. We assessed normalized current signals, average base quality scores, and dwell times at CpG motifs (Figures 1F,G). Across all three metrics, clear differences were observed between methylated and unmethylated samples—most notably in current signal profiles. Despite these distinctions, the noisy nature of the data again hindered the effectiveness of conventional statistical techniques in drawing reliable boundaries between methylated and unmethylated states. These limitations strongly motivated the adoption of deep learning-based models for simultaneous base calling and methylation detection.

We introduce Bream, an open-source deep learning framework that enables simultaneous base calling and DNA methylation detection. As depicted in Figure 2A, Bream integrates a comprehensive workflow, processing raw signals to output both sequence data and methylation status. The framework’s core consists of two primary components: a base calling model based on a convolutional neural network (CNN) combined with a reverse long short-term memory (LSTM) network (Figure 2B), and a methylation detection model leveraging a bidirectional LSTM (BiLSTM) with an attention mechanism (Bahdanau et al., 2016) (Figure 2C). Bream achieves high accuracy in decoding base sequences and identifying CpG site methylation directly from nanopore raw signals.

This design ensures a seamless pipeline from input to output, emphasizing clarity and efficiency in handling complex genomic data. By streamlining these processes, Bream advances the capability for integrated analysis in bioinformatics applications.

We systematically evaluated the performance of the Bream framework using datasets from multiple species, including A. thaliana, O. sativa, and D. melanogaster with defined methylation states (methylated and unmethylated, treated with methyltransferase (Abdelraheem et al., 2022)). The Bream model, trained on these datasets, achieved efficient base calling across all test datasets (Supplementary Table 5). We then analyzed the quality of the sequencing data and the accuracy of the base calling.

In assessing dataset quality, we calculated the average base quality and evaluated read and base pass rates within the Fastq files of each dataset. Results, depicted in Figures 3A,B, highlight A. thaliana as the highest quality dataset, with a read pass rate of 63.67% and a base pass rate of 84.17%. In contrast, O. sativa demonstrated inferior quality, evidenced by a read pass rate of 49.79% and a base pass rate of 71.21%. Figure 3C further illustrates the quality distribution of A. thaliana reads, revealing a notable proportion of low-quality reads characterized by an average quality score ranging from 0 to 2. To evaluate base calling error rates, we employed minimap2 (Li, 2022) to align reads (with an average base quality ≥10) to reference genomes, thereby extracting alignment accuracy (identity) from the CIGAR information. As depicted in Figures 3A,D thaliana achieved the highest alignment accuracy at 91.83%, while the YF6418 sample recorded the lowest alignment accuracy at 89.38%. The error, deletion, insertion, and mismatch rates for A. thaliana are detailed in Figure 3E, with comprehensive results from the full dataset provided in Supplementary Table 1, indicating a nanopore sequencing error rate ranging from 8% to 10%.

Finally, the accuracy of the assembled sequences derived from base calling was assessed using the software Merqury (Rhie et al., 2020), which evaluates phred base quality values. As illustrated in Figure 3F, there is a correlation between assembly accuracy and alignment accuracy across the various datasets.

We trained methylation calling model on D. melanogaster dataset. Evaluation on the D. melanogaster dataset included calculating key metrics such as ROC and PR curves (Figures 4A,B), yielding areas under the curve of 98.35% and 98.79%, respectively. Supplementary Table 3 details the model’s precision, recall, and F1 scores at 96.79%, 93.88%, and 95.31%, respectively.

For further validation of real-world data, we analyzed DNA samples from A. thaliana, and O. sativa leaves using bisulfite sequencing data and nanopore sequencing data. CpG sites with 100% or 0% methylation frequency in the bisulfite data were used as ground truth labels for binary classification. Performance metrics showed that for A. thaliana, accuracy, recall, and F1 scores were 94.81%, 92.74%, and 94.70% (Supplementary Table 3), respectively. For Oryza, these metrics were 97.38%, 97.33%, and 97.38% (Supplementary Table 3). These results demonstrate that the Bream methylation model has a strong capability for single-molecule CpG methylation detection.

At the genomic level, Bream’s methylation predictions were compared with whole genome BS-seq methylation frequencies by analyzing their correlation across sequencing depths. As shown in Figures 4C,D, the quantitative correlation improved with increasing sequencing depth but plateaued beyond 15× coverage. Heatmaps of A. thaliana and O. sativa datasets (Figures 4E,F) revealed high consistency between Bream predictions and BS-seq data. In the A. thaliana dataset (N = 5,488,306), Pearson’s correlation coefficient was 0.974, while for O. sativa (N = 29,503,269), it was 0.977. These results highlight the substantial agreement between Bream and BS-seq in methylation detection while demonstrating Bream’s ability to cover more CpG sites, underscoring its significant advantages in methylation analysis.

The accuracy of base calling in nanopore sequencing is highly dependent on signal quality. The technology is susceptible to high levels of noise, and certain sequencing regions may exhibit poor signal characteristics. These high-noise regions, characterized by elevated signal fluctuations, adversely affect read accuracy and the reliability of methylation detection. Specifically, in base calling, high-error regions indicate inherently unreliable reads, which directly compromise downstream analyses such as genome assembly and variant detection. For methylation detection, which relies on discerning subtle signal patterns to identify modified bases, high-error regions lead to misalignment of methylation sites. This introduces inaccuracies in training data, ultimately impairing prediction performance. To mitigate these issues, computational strategies such as error correction or filtering of low-quality regions have been developed. For example, NECAT (Chen et al., 2021) employs an adaptive algorithm to correct errors in de novo assembly, while DeepSignal-plant (Ni et al., 2021) refines modification calls by aligning sequences to a reference genome to rectify erroneous k-mers. Nevertheless, the most critical factors for improvement remain continual advancements in sequencing chemistry, flow cell design, and the precision of base calling and methylation detection algorithms.

To evaluate the data quality from our Bream output, we filtered out reads with an average sequencing quality score below 10 across four datasets. We then examined the distribution of average read errors after filtering (Figure 5A). The results show that base calling errors are generally distributed between 0.05 and 0.15, with noticeable variation across different samples. The error distribution for A. thaliana is skewed to the left, indicating lower error rates, whereas the error rates for the synthetic D. melanogaster datasets (YF6418 and YF6419) are similar. Notably, YF6418, treated with methyltransferase, shows a slight right skew, suggesting that methylation may affect sequencing accuracy.

Furthermore, we explored the relationship between Bream-reported base quality scores and the read accuracy estimated using consensus sequences of corresponding Pacbio HiFi data. Using a 500 bp window, we compared average sequencing quality scores with error rates after aligning the reads to the reference genome (Figure 5B). The result illustrates that as the quality score increases, the error rate decreases significantly. When the average quality exceeds 25, error rates approach zero, indicating a potential correlation between higher sequencing quality and fewer errors.

To evaluate the impact of regions with high signal noise on methylation detection, we analyzed the methylation prediction accuracy under different error rates in the D. melanogaster in vitro synthetic datasets (Figure 5C). The results reveal that low error rates (<0.1) lead to high prediction accuracy, while error rates above 0.25 cause a significant drop in accuracy to around 70%. This underscores the importance of maintaining low error rates to ensure reliable methylation detection. Excluding regions with high error rates may improve prediction accuracy by minimizing the impact of poor-quality data.

Finally, we plotted ROC curves and calculated AUC values for methylation prediction at different error rate ranges (Figure 5D). At low error rates (e.g., 0–0.05 or 0.05–0.1), the AUC values approach 1, indicating excellent prediction accuracy. As error rates increase, the ROC curves shift downward, and AUC values decrease. Notably, at error rates ≥0.3, the AUC drops to 0.79, further emphasizing that higher error rates significantly reduce the accuracy of methylation prediction.

Recombinant plasmids containing sequences of either Pifi-like helicase or mutated pore proteins (amino acid sequences detailed in Supplementary Table 7) were transformed into BL21 (DE3) competent cells using a heat-shock protocol. Following transformation, the cells were plated on solid LB agar containing ampicillin and incubated overnight at 37 °C. Single colonies were selected and cultured in liquid LB medium supplemented with ampicillin, shaken at 200–220 rpm at 37 °C. Optical density (OD600) was periodically measured to monitor cell growth.

Protein expression was induced upon reaching specific optical densities: OD600 of 0.6–0.8 for helicase cultures and OD600 of 2.0–2.2 for pore protein cultures. The cultures were cooled to 16 °C–18 °C, and IPTG was added at final concentrations of 1 mM for helicase or 0.015 mM for pore proteins. Induction was maintained for 12–24 h.

Post-induction, cells were harvested via centrifugation, lysed by high-pressure homogenization, and the target proteins were purified. Pifi-like helicase was purified using FPLC, while pore proteins were isolated through Ni-NTA affinity chromatography. Eluted protein samples were collected for further analysis.

We generated initial DNA sequences using the base calling model provided by Qitan Technology. These sequences were aligned to the original electrical signals, and corrections were made using the Minimap2 tool to improve accuracy. During preprocessing, we first removed invalid signals representing the stage before DNA entered the nanopore using signal detection methods. Next, the signals were normalized using the Median Absolute Deviation (MAD) (Ni et al., 2019) (Bai et al., 2024) approach. The data was then divided into segments of 6,000 base pairs (bp) with a 500 bp overlap. Chunks with coverage less than or equal to 95% were filtered out to ensure high-quality input data for base calling model training.

The Bream framework comprises six computational modules (Figure 2A): signal preprocessing, basecalling, CTC (Graves et al., 2023) decoding, consensus assembly, FASTQ quality value conversion and writing, and methylation detection. To ensure consistent basecalling, the preprocessing module divides raw signals of varying lengths into overlapping segments, each containing 6,000 sample points. These segments are then processed by the basecalling module to identify nucleotide signals. The CTC decoding module interprets these signals to reconstruct the base sequence. Subsequently, the consensus assembly module integrates overlapping segments from the same read to reconstruct the full sequence along with its corresponding quality scores. The writing module formats this information into a standard FASTQ file. Simultaneously, the methylation detection module analyzes CpG motifs within the reconstructed sequences and reports their methylation status.

The base calling module of Bream (Figure 2B) employs a five-layer convolutional neural network (CNN) followed by a five-layer bidirectional long short-term memory (BiLSTM) network. Input signal segments are processed through convolutional layers and batch normalization to extract local signal features and standardize variations. The resulting features are then passed to the bidirectional LSTM layers, which capture contextual information from both forward and reverse directions. A linear layer subsequently converts these features into probabilities for each of the four nucleotide bases (A, C, G, T) or a blank symbol at every position. Finally, a connectionist temporal classification (CTC) decoder translates the probability sequence into the final base sequence. The CTC approach enables unsupervised alignment by allowing variable-length output sequences without strict signal-to-base correspondence, while probabilistically merging paths to resolve input-output length discrepancies.

The training dataset comprised synthetic data from A. thaliana, O. sativa, and D. melanogaster (datasets YF6418 and YF6419), with a sampling ratio of 4:2:2:2. Approximately 5 million signal segments were extracted for training. Training was conducted on a server equipped with an AMD 9654 CPU and two NVIDIA L40s GPUs. A batch size of 256 was used with the AdamW (Loshchilov and Hutter, 2019) optimizer, an initial learning rate of 0.001, and a linear warm-up phase followed by cosine decay scheduling. To accelerate training, PyTorch’s automatic mixed precision was employed. The entire training process spanned approximately 80 epochs over 3 days. During inference, the model was converted to half-precision for GPU execution, significantly improving throughput (measured in bases called per second).

The methylation detection module (Figure 2C) employs a bidirectional long short-term memory network augmented with an attention mechanism. This module integrates two types of features: (1) sequence-based features, represented as one-hot encodings of the 21-base window centered on each CpG site, and (2) normalized raw signal features derived from the corresponding genomic region. These features are processed independently through separate BiLSTM layers to capture methylation-related patterns, after which their outputs are concatenated. An attention layer dynamically weights the contributions of these features, followed by additional BiLSTM and linear layers that compute probabilities for unmethylated (P _um) and methylated (P _m) states at each CpG site.

The Bream methylation model was trained with fully methylated (CpG methyltransferase-treated) and fully-unmethylated (PCR-treated) D. melanogaster datasets. These annotated datasets were balanced with an equal ratio of positive and negative samples (1:1), randomly partitioned into 90% for training and validation, and 10% for testing model performance. To train the model, methylation-related features were extracted using basecall sequences generated by the Bream model and alignment data. The model was trained on approximately 100 GB of synthetic data from D. melanogaster, optimized with the Adam optimizer (Kingma and Ba, 2017) at a learning rate of 0.001 and a linear decay factor of 0.4.

To assess dataset quality and model performance, we systematically evaluated base calling results. Only reads with an average base quality value (Q-score) ≥10 were retained for downstream analysis, a threshold selected based on previous studies (Zhang et al., 2021) and the overall performance of Bream on the Qitan platform. Raw data quality was evaluated by calculating the average base quality, read pass rate, and base pass rate from the FASTQ files of each dataset. To assess base calling error rates, filtered reads were aligned to the reference genome using minimap2. Alignment accuracy (identity) along with error, deletion, insertion, and mismatch rates were derived from the CIGAR information. Finally, the accuracy of the assembled sequences generated through base calling was evaluated using Merqury.

For methylation data, in addition to presenting classification accuracy using ROC curves on a methyltransferase-treated Drosophila dataset, we calculated the Pearson correlation coefficient between per-site methylation frequencies derived from our method and those obtained from second-generation bisulfite sequencing in Arabidopsis and rice. The coefficient is defined as: ρ X , Y = C o v X , Y σ X σ Y

Where C o v denotes covariance, and σ X and σ Y represent the standard deviations of data X and Y, which correspond to per-site methylation frequencies calculated by our method and bisulfite sequencing, respectively.