Accurate detection of tandem repeats from error-prone sequences with EquiRep

Evaluation with simulated random sequences

The simulated random sequences are generated as follows: (1) generate a random string U constituting nucleotides (A,T,G,C) of length five, 10, 50, 100, 200, 500, and 1000 that serves as the ground-truth repeat unit; (2) concatenate multiple copies of the unit U to generate a longer sequence, with the frequency (number of copies) of the unit being three, five, 10, and 20; (3) introduce random errors—insertions, deletions, and substitutions at equal probabilities—at the rates of 10%, 15% and 20% into the concatenated string to simulate real-world sequencing errors and mutations; and (4) insert random strings, matching the length of the concatenated string (i.e., the repeat region), at both sides of the concatenated string.

For each of the settings (the combination of unit length, frequency of units, and error rate), we randomly and independently generated 50 sequences. We evaluated the methods’ predictions as follows. Let T be a ground-truth repeat unit, and let P be a prediction. We computed a rotation-aware edit distance between P and T. Specifically, because P may be a rotation of the T, we calculated the edit distance between T and all possible rotations of P and take the minimum value, defined as the rotation-aware edit distance. For each setting, we analyzed the 50 instances and report the following three metrics. First, we measured accuracy as the number of instances (out of 50) in which the method predicts the exact ground-truth unit (i.e., rotation-aware edit distance is zero). Second, we evaluated the proportion of close predictions, defined as cases in which the rotation-aware edit distance is <10% of the true unit length. Third, we reported the average of the normalized rotation-aware edit distance (distance divided by the unit length) across all 50 instances.

Figure 1, A–G, compares the accuracy on simulated data at a 10% error rate for various lengths and copy numbers. EquiRep consistently predicts a comparable or greater number of correct instances compared with other methods. The methods with a performance closest to that of EquiRep appear to be mTR and TRF; however, both struggle to maintain accuracy with large unit lengths. The accuracy of EquiRep is significantly higher than any of the other methods for unit lengths of 500 and 1000 bp, which demonstrates the ability of our tool to predict longer tandem repeats. Figure 2, A–G, compares the ratio of close predictions on simulated data at a 10% error rate. The ratio for EquiRep is high regardless of the copy number, and the trend tends to be consistent over the different unit lengths, unlike other methods. Figure 3, A–G, compares the averaged normalized rotation-aware edit distance. Observe that EquiRep consistently achieves the lowest distance, indicating that even when its predictions are incorrect, they remain the closest to the true sequence.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

Figure 1.

Comparison of accuracy on simulated data at a 10% error rate, for (A–G) unit lengths 5, 10, 50, 100, 200, 500, and 1000, respectively.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

Figure 2.

Comparison of proportion of close predictions (rotation-aware edits <10% of the unit length) on simulated data at a 10% error rate, for (A–G) unit lengths 5, 10, 50, 100, 200, 500, and 1000, respectively.

View larger version:

• In this window
• In a new window

• Download as PowerPoint Slide

Figure 3.

Comparison of average normalized rotation-aware edit distance on simulated data at a 10% error rate, for (A–G) unit lengths 5, 10, 50, 100, 200, 500, and 1000, respectively.

To better illustrate the distributions of the normalized rotation-aware edit distances between the predicted unit and the ground-truth, we show the fine-grained plots for all simulated settings on data with a 10% error rate, available in Supplemental Figures S5–S31.

A comparison with another approach, dot2dot (Genovese et al. 2019), on simulated data with 10% error rate, is available at Supplemental Figure S32, A–G, through Supplemental Figure S34, A–G. EquiRep outperforms dot2dot drastically on all settings.

Results for 15% and 20% error rates are available in Supplemental Figure S35, A–G, through Supplemental Figure S37, A–G, and in Supplemental Figure S38, A–G, through Supplemental Figure S40, A–G, respectively. For a higher error rate, TRF, mreps, and TideHunter see a sharp decline in accuracy as the unit length exceeds 10 bp. Conversely, mTR's ability to handle long, noisy reads allows it to achieve an accuracy close to that of EquiRep; however, its performance drops when the unit length reaches ≥500 bp. At such a long unit and with high sequencing errors, all methods struggle to accurately predict tandem repeats, particularly when the copy number is low. Overall, EquiRep outperforms other tools on the three metrics across different simulations.

Evaluation with data simulated with PBSIM2

To better mimic real long reads, we evaluated our method using data simulated by PBSIM2 (Ono et al. 2021). To simulate, we first generated sequences containing repeats positioned in the middle with random sequences flanking both ends. The repeat configurations were consistent with those described in the section “Evaluation with simulated random sequences,” including repeat units of lengths five, 10, 50, 100, 200, 500, and 1000, with each unit repeated three, five, 10, or 20 times. The following command (pbsim ‐‐depth 1 ‐‐hmm_model PC64.model ‐‐accuracy-mean 0.90) is subsequently used to simulate long reads using PBSIM2. Results were compared against the same set of alternative methods, detailed in Supplemental Figures S41, A–G, through Supplemental Figure S43, A–G. EquiRep consistently outperformed the competing methods in nearly all scenarios, highlighting its effectiveness on more realistic simulated reads.

Evaluation using simulated sequences with recurring k-mers in a unit

Genomic sequences are not purely random, often containing recurring substrings. We compared different methods on this scenario with simulations in which the repeat unit itself contains recurring structures. In this setting, predicting the correct repeat sequence is challenging as methods may encounter difficulties in distinguishing between such recurring k-mers in a single unit and identical k-mers across multiple units.

We used this approach to simulate the above sequences. First, for a given unit length l ∈ {50, 200, 500}, we generated a random k-mer of length k ∈ {5, 10, 20}, respectively. Second, we constructed the repeat unit by concatenating the random k-mer two or three times. After these concatenations, any remaining positions within the unit (i.e., l − 2k for two concatenations and l − 3k for three concatenations) will be filled with random nucleotides. Third, we concatenated multiple copies of the repeat unit to generate a longer sequence, with the frequency of units being three, five, 10, and 20. Fourth, we introduced random errors at rates of 10% and 20%. Fifth, at the end, we inserted random strings, matching the length of the concatenated string at both ends.

The same evaluation metrics for the previous simulations are also used here. Supplemental Figure S44, A–C, indicates accuracy (the ratio of fully correct instances) of EquiRep exceeding or equaling the other methods when the simulations have two copies of a k-mer within the unit at a 10% error rate. Supplemental Figure S45, A–C, shows that almost for all instances the edits predicted by our method are <10% of the unit length. Again, EquiRep achieves the lowest averaged distance as illustrated in Supplemental Figure S46, A–C. Supplemental Figures S47, A–C, through S49, A–C, demonstrate the results for data with a 20% error rate. There is a drastic decline in accuracy for all methods except mTR and EquiRep.

Nearly all repeat units generated by EquiRep have edits <10% of the unit length for copy numbers above 10, which highlights the reliability of our predictions, especially in challenging erroneous settings.

We also tested all methods on another set of data with three copies of repeating k-mers within the repeat unit, as shown in Supplemental Figure S50, A–C, through Supplemental Figure S52, A–C (for an error rate of 10%) and in Supplemental Figure S53, A–C, through Supplemental Figure S55, A–C (for an error rate of 20%). EquiRep is able to make better or similar predictions in all cases, indicating that its algorithm is the least affected by the presence of embedding k-mers within repeat units.

Evaluation using human satellite DNA data

We then tested all methods on reconstructing repeat units for satellite DNA in human Chromosome 5 (Paar et al. 2007). This known satellite DNA consists of 13 units (i.e., 13 monomers), each of which is ∼171 bp in size. To construct the input sequence for methods to predict, we concatenated the 13 monomers into a string denoted as (x). To create more testing instances, we introduced flanking regions on both sides of the concatenation denoted as (axa) and introduced errors of 1%, 5%, and 10% to (x) and (axa). To evaluate the predicted unit by different methods, we calculated the normalized rotation-aware edit distance between the predicted unit with each of the 13 known monomers and reported the averaged distance.

Table 1 shows the results. EquiRep consistently maintains a lower normalized distance, outperforming or matching all other tools. The values for EquiRep are similar to those of mTR when the input sequences have flanking regions at either end (axa), but our method is ∼87% better than mTR when just the repeat region is provided (x). Although TideHunter and TRF exhibit accuracy levels similar to ours, they fall short at higher error rates, for which EquiRep excels, with an 87% improvement.

Evaluation using Caenorhabditis elegans centromere ONT data

We adopted a data set reported by Yoshimura et al. (2019), who studied the assembly of the Caenorhabditis elegans genome using Nanopore long-read data. We collected the raw long reads that are aligned to the centromere (listed in their Supplemental Fig. S4). Each of the long reads may contain more than one repeating region. Because our current method does not support detecting multiple repeating regions in a single input sequence, we manually extracted the rough region with repeats. Specifically, we first generated a dot plot for each long read, observed the repeating regions, and then manually cut out these regions and piped them to each of the methods. The ground-truth sequence of the unit is available, which is obtained by curating from PacBio HIFI data sets.

Table 2 presents the normalized rotation-aware edit distance between the predicted units and the ground truth. We report the average value across all cases. EquiRep achieves the second-best performance. For each method, we also report the number of cases in which the normalized rotation-aware edit distance is below 0.2, indicating high-quality predictions. EquiRep performs well in seven out of 13 cases, whereas the top-performing methods, mTR and TRF, achieve good predictions in eight cases.

Evaluation with RCA data

The set of real data is a RCA-based ONT sequencing protocol from isocirc (Xin et al. 2021) that has been used to detect a catalog of full-length circular RNAs from 12 human tissues. We considered a subset of 101 sequences from prostate tissue long-read ONT data (obtained from the NCBI Gene Expression Omnibus [GEO; https://www.ncbi.nlm.nih.gov/geo/] under accession number GSE141693) for analysis. It is difficult to evaluate the repeats from the RCA-based long-read data owing to lack of reliable ground truth, so we evaluated these data in two different ways. First, we used a dot plot analysis. Dot plots have served as a common approach for visualizing and identifying the structural patterns of sequences such as repeats. We first aligned the input sequence to itself with LASTZ (Harris 2007) using specific parameters designed for generating dot plots. The alignment program generates a dot file that can be converted to an image file for visualization using a simple R (R Core Team 2021) script. The dot file can be used to estimate the repeat unit length (but not sequence of the unit). We treated this estimate as a benchmark for comparing the predictions of EquiRep and other tools. We reported the number of predictions that fall within 5%, 20%, 50%, and 80% error range of the true length. For the second approach, we first concatenated copies of the unit predicted to get a string A, which is longer than the input sequence. Then we got the “semiedit-distance,” which is the smallest edit distance between any substring of A and the input sequence. The idea behind this metric is that if the prediction is accurate, then the multiple concatenation of it should match the input sequence very well. We recorded the smallest edit distance and report the number of instances on which a method has a ratio (semiedit-distance)/(input sequence length) less than or equal to 0.1, 0.2, 0.3, 0.5, and 0.8.

Table 3 compares different methods in terms of the predicted repeat unit length, and Table 4 compares the normalized semiedit-distance. In both metrics, EquiRep demonstrates high accuracy, consistently outperforming mTR, TRF, and mreps. The results are also comparable to TideHunter, which is specifically optimized for RCA-based analysis. Given that the exact repeat sequences for this data set are not available, similar metric values in the table can be interpreted as comparable accuracy. It should be noted that although TideHunter excels on RCA data, its accuracy diminishes on shorter unit repeats, as indicated by the simulation results. This highlights that EquiRep is adaptable to a broad range of complex sequences and is versatile for various applications.

In the above analysis of the RCA data sets, we observed that many repeat units exceed 1000 bp in length. This is consistent with the fact that many expressed circular RNAs are themselves >1000 bp. These observations also support the use of longer unit lengths (e.g., 500 bp and 1000 bp) in our simulated experiments (see section “Evaluation with simulated random sequences”).

Analysis of sensitivity of EquiRep to parameters

We conducted experiments to analyze the sensitivity of EquiRep to its three key parameters: (1) the score threshold (default: 25) used to identify significant paths from the initial matrix D; we tested the alternative values zero, 10, and 50; (2) the window size (default: 7) used for identifying local maxima in initial matrix D; we tested two other choices, five and nine; and (3) the number of iterations (default: 5) of iterative matrix refinement; we tested two other values, one and 10. To assess the effect of a choice of a parameter, we made it the only change to the default setting of EquiRep and then compared the variant with the default EquiRep. The same simulated data, used in the section “Evaluation with simulated random sequences,” with a 10% error rate was used here to obtain the results. We also used the same three metrics in the evaluation.

The results corresponding to the three parameters were given, respectively, in Supplemental Figure S56, A–G, through Supplemental Figure S58, A–G, Supplemental Figure S59, A–G, through Supplemental Figure S61, A–G; and Supplemental Figure S62, A–G, through Supplemental Figure S64, A–G. We can conclude that EquiRep is not sensitive to any of them, justifying its default choices.

Comparison of running time

Supplemental Table S1 presents the runtime of all methods on the simulated data from the section “Evaluation with simulated random sequences” with a 10% error rate. On average, mTR had the longest runtime, followed by EquiRep. TRF, mreps, and TideHunter were significantly faster. As noted in the Discussion, several modules in EquiRep are parallelizable, and we are optimistic about further improving its computational efficiency.

EquiRep is well-suited for processing a large number of error-prone long reads on multicore servers, as it operates on individual reads, allowing efficient batch processing that fully utilizes available cores. It is also likely that, in large-scale long-read data set, the majority of the long reads do not contain repeating regions. Fast filtering strategies, such as the seed-chaining procedure used in step 1 of EquiRep, can quickly discard such reads, leaving only a small subset that requires full processing by the complete EquiRep algorithm.