Error assembly detection and repair method and system based on double-end sequencing data
Through error assembly detection and repair methods based on double-ended sequencing data, the characteristic signals of candidate error regions are extracted and clustered to analyze the characteristics of candidate error regions, and the assembly errors are identified and repaired, which solves the problem of frequent assembly errors in the prior art, and improves the accuracy and completeness of genome assembly.
Patent Information
- Application Number
- CN202510517973.0
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-04-24
- Publication Date
- 2025-05-27
- Estimated Expiration
- 2045-04-24
AI Technical Summary
Existing genome assembly algorithms are prone to errors during assembly, resulting in structural errors such as insertion/deletion errors and incorrect connection errors, affecting the accuracy of the genome and the quality of downstream research.
Using error assembly detection and repair methods based on double-ended sequencing data, by comparing the assembly results with the double-ended sequencing data, a variety of characteristic signals of candidate error areas are extracted, such as coverage abnormal signals, insertion or missing abnormal signals, and cluster analysis is carried out to identify the error assembly area and repair it according to the error type.
Effectively identify and repair assembly errors, improve the accuracy and completeness of genome assembly, reduce false positive rates, be able to mark the specific location of the error, and have the function of correcting assembly errors.
Abstract
Description
Technical Field
[0001] The present invention belongs to the technical field of sequence analysis, and particularly relates to a method and system for detecting and repairing incorrect assemblies based on paired-end sequencing data. Background Art
[0002] With the rapid development of genome sequencing technology, various assembly algorithms have emerged one after another, greatly improving the ability to construct genomic sequences from sequencing data. Complete and accurate genome assembly is crucial for analyzing the genetic information of organisms, studying evolutionary mechanisms, and identifying gene variations related to diseases (including single nucleotide polymorphisms, insertions / deletions, and structural variations). This information has important application values in fields such as disease mechanism research, precision medicine, and personalized treatment.
[0003] However, achieving high-quality assembly is challenging. Existing genome assembly algorithms still have certain limitations, which may lead to errors during the assembly process. Incorrect assemblies usually manifest as structural errors, such as insertion / deletion (indel) errors and misjoin errors. These errors not only affect the accuracy of the genome but also have a negative impact on comparative genomics, genetic variation analysis, and downstream bioinformatics research. To solve the problem of incorrect assemblies, various error detection methods have been proposed. These methods are mainly divided into two categories: reference genome-based methods and read alignment-based methods.
[0004] Reference genome-based methods detect errors by aligning the assembly results with the reference genome. For example, QUAST can identify errors such as translocations and inversions, but it is difficult to distinguish between assembly errors and true variations generated by individuals. In addition, due to the lack of accurate reference genomes for many organisms in the real environment, the applicability of this type of method is limited. Read alignment-based methods, on the other hand, use the sequencing data itself to evaluate the assembly quality, and identify repetitive regions, incorrect assemblies, and indel errors by analyzing features such as changes in coverage depth and abnormal insert fragment lengths. However, existing read alignment-based methods have a high false positive rate, and most tools cannot mark the specific positions where errors occur, let alone have the function of correcting assembly errors, which affects their application in actual scenarios. Summary of the Invention
[0005] Current sequencing technologies and genome assembly tools still face many challenges in the process of complete chromosome reconstruction. These challenges mainly include sequencing errors, limitations in sequencing coverage, and a large number of repetitive sequences present in the genome. Although second-generation sequencing (NGS) technologies provide short-read data with high coverage and high accuracy, due to factors such as polyploid structure, genetic diversity, heterozygosity, and repetitive sequences, errors still commonly exist in genome assembly. Therefore, the present invention provides a method and system for detecting and repairing misassemblies based on paired-end sequencing data.
[0006] The technical solution of the present invention is as follows: The present invention provides a method for detecting and repairing misassemblies based on paired-end sequencing data, including the following steps: S1: Based on the assembly result, with a preset suspicious site as the core, extend a first preset length to both sides of the genomic sequence, and use it as a candidate error region; Compare the paired-end sequencing data with the assembly result. Based on the comparison information, extract various characteristic signals of the candidate error region. The various characteristic signals are coverage anomaly signals, insertion or deletion anomaly signals, paired direction anomaly signals, insertion size anomaly signals, read pair anomaly signals, and splicing anomaly signals; S2: For each candidate error region, based on the various characteristic signals of the candidate error region extracted, after clustering, obtain the clustering result of each characteristic signal of the corresponding candidate error region; if the clustering result of one characteristic signal is abnormal, then this candidate error region is a misassembly region; S3: Based on the misassembly region, determine the type of assembly error. If the type is insertion or deletion, remove the inserted sequence or the inserted and deleted fragment; If the type is incorrect connection, generate a consensus sequence through multiple sequence alignment, compare the consensus sequence with the assembly result, and perform repair.
[0007] In S1, the specific method for extracting the coverage anomaly signal of the candidate error region is as follows: Taking the candidate error region as the benchmark, extend a second preset length to both sides of it, and use it as an extended region; Calculate the average coverage of all sites in the extended region, traverse and detect the coverage of each site in the candidate error region. If the coverage of a certain site is not less than the first preset multiple of the average coverage or not greater than the second preset multiple of the average coverage, then this site is marked as a coverage anomaly site; Count the total number of coverage anomaly sites in the candidate error region, calculate the ratio of the total number to the total length of the candidate error region, and use this ratio as the coverage anomaly signal.
[0008] In S1, the specific method for extracting the insertion or deletion anomaly signal of the candidate error region is as follows: Statistically count the number of reads with insertions or deletions at each locus in the candidate error region. If the coverage of a certain locus is not less than the preset value, calculate the ratio of the number of reads to the coverage of this locus, and use this ratio as the insertion or deletion signal of this locus; Compare the insertion or deletion signals of all loci in the candidate error region, and use the maximum ratio as the insertion or deletion abnormal signal of this candidate error region.
[0009] In step S1, extract the abnormal signal of the paired direction in the candidate error region, specifically: Detect reads with mismatched directions at each locus in the candidate error region. If reads with mismatched directions are found at a certain locus, mark the local area of this locus as the region with an incorrect paired direction, and statistically count the number of reads with mismatched directions and the total number of reads in the region with an incorrect paired direction. Use the ratio of the number of reads to the total number as the abnormal signal of the paired direction in this region with an incorrect paired direction; Compare the magnitudes of the abnormal signals of the paired directions in all regions with incorrect paired directions within the candidate error region, and use the abnormal signal of the paired direction in the region with the maximum value as the abnormal signal of the paired direction in this candidate error region.
[0010] In step S1, extract the abnormal signal of the insertion size in the candidate error region, specifically: Detect reads with abnormal insertion sizes at each locus in the candidate error region. If reads with abnormal insertion sizes are found at a certain locus, mark the local area of this locus as the region with an abnormal insertion size, and statistically count the number of reads with abnormal insertion sizes and the total number of reads in the region with an abnormal insertion size. Use the ratio of the number of reads to the total number as the abnormal signal of the insertion size in this region with an abnormal insertion size; Compare the magnitudes of the abnormal signals of the insertion sizes in all regions with abnormal insertion sizes within the candidate error region, and use the abnormal signal of the insertion size in the region with the maximum value as the abnormal signal of the insertion size in this candidate error region.
[0011] In step S1, extract the abnormal signal of the read pair in the candidate error region, specifically: Detect unpaired reads at each locus in the candidate error region. If unpaired reads that cannot be aligned to the same assembly result are found at a certain locus, mark the local area of this locus as the region of unpaired reads, and statistically count the number of reads with abnormal pairings and the total number of reads in the region of unpaired reads. Use the ratio of the number of reads to the total number as the abnormal signal of the read pair in the region of unpaired reads; Compare the magnitudes of the abnormal signals of the read pairs in all regions of unpaired reads within the candidate error region, and use the abnormal signal of the read pair in the region with the maximum value as the abnormal signal of the read pair in this candidate error region.
[0012] S1 extracts the splicing abnormal signal of the candidate error region, specifically: Count the number of reads with splicing at each site in the candidate error region one by one. If the coverage of a certain site is not less than the preset value, calculate the ratio of the number of reads to the coverage of this site, and take this ratio as the splicing abnormal signal of this site; Compare the splicing abnormal signals of all sites in the candidate error region, and take the maximum ratio as the splicing abnormal signal of this candidate error region.
[0013] Before S2 clustering, it further includes: Randomly select several regions from the assembly results longer than the preset sequence length, and jointly perform abnormal signal clustering analysis with the candidate error regions.
[0014] In S3, if the type is incorrect connection, generate a consensus sequence through multiple sequence alignment, and compare the consensus sequence with the assembly result for repair, specifically: Based on the assembly result of the incorrect connection type, analyze the alignment characteristic signals and mark the breakpoint region; Based on the paired-end sequencing data, extract the paired-end reads that completely cover the breakpoint region, and after filtering, obtain the local read sequence; The local read sequence generates a consensus sequence through multiple sequence alignment. Compare the consensus sequence with the assembly result. If the successfully aligned length of the repaired sequence accounts for more than 90% of the breakpoint region, replace the incorrect connection region in the assembly result with the consensus sequence.
[0015] The present invention also provides an incorrect assembly detection and repair system based on paired-end sequencing data, including: Feature signal extraction module: Based on the assembly result, with the preset suspicious site as the core, extend a first preset length to both sides of the genomic sequence, and use it as the candidate error region; Compare the paired-end sequencing data with the assembly result, and based on the comparison information, extract various feature signals of the candidate error region. The various feature signals are coverage abnormal signal, insertion or deletion abnormal signal, paired direction abnormal signal, insertion size abnormal signal, read pair abnormal signal, splicing abnormal signal; Incorrect assembly region determination module: For each candidate error region, after clustering based on the various feature signals of the extracted candidate error region, obtain the clustering result of each feature signal of the corresponding candidate error region; if there is an abnormal clustering result of one feature signal, then this candidate error region is an incorrect assembly region; Repair module: Based on the incorrect assembly region, determine the type of assembly error. If the type is insertion or deletion, remove the inserted sequence or the inserted and deleted fragment; If the type is incorrect connection, a consensus sequence is generated through multiple sequence alignment, and the consensus sequence is aligned with the assembly result for repair.
[0016] Beneficial effects The present invention utilizes the alignment information of paired-end sequencing data with the assembly result, extracts various abnormal signals, and performs clustering analysis, which can effectively identify true misassemblies; according to the type of misassembly, combined with the sequencing data information, the assembly result is corrected, thereby optimizing the assembly result and improving the overall accuracy and integrity of genome assembly. Specific implementation manners
[0017] The following embodiments are intended to illustrate the present invention rather than further limit the present invention.
[0018] The present invention provides a method for detecting and repairing misassemblies based on paired-end sequencing data, including the following steps: S1: Based on the assembly result, with a preset suspicious site as the core, extend a first preset length to both sides of the genomic sequence to obtain a candidate error region; The paired-end sequencing data is aligned with the assembly result, and based on the alignment information, various characteristic signals of the candidate error region are extracted. The various characteristic signals are coverage abnormal signals, insertion or deletion abnormal signals, paired direction abnormal signals, insertion size abnormal signals, read pair abnormal signals, and splicing abnormal signals; S2: For each candidate error region, based on the various characteristic signals of the candidate error region extracted, after clustering, the clustering result of each characteristic signal of the corresponding candidate error region is obtained; if the clustering result of one characteristic signal is abnormal, then the candidate error region is a misassembly region; S3: Based on the misassembly region, the type of assembly error is determined. If the type is insertion or deletion, the inserted sequence or the inserted and deleted fragment is removed; If the type is incorrect connection, a consensus sequence is generated through multiple sequence alignment, and the consensus sequence is aligned with the assembly result for repair.
[0019] The present invention utilizes the alignment information of paired-end sequencing data with the assembly result (contig), extracts various abnormal signals, and performs clustering analysis, which can effectively identify true misassemblies; according to the type of misassembly, combined with the sequencing data information, the assembly result is corrected, thereby optimizing the assembly result and improving the overall accuracy and integrity of genome assembly.
[0020] The following describes each step of the operation in detail.
[0021] S1: Based on the assembly result, with a preset suspicious site as the core, extend a first preset length to both sides of the genomic sequence, and use it as the candidate error region (for example, the total length of the candidate error region is 200bp); Align the paired-end sequencing data with the assembly result. Based on the alignment information, extract multiple characteristic signals of the candidate error region. The multiple characteristic signals are coverage anomaly signals, insertion or deletion anomaly signals, paired direction anomaly signals, insertion size anomaly signals, read pair anomaly signals, and splicing anomaly signals.
[0022] Considering that: (1) Correctly assembled sequences usually exhibit uniform coverage. Deviations from the expected coverage (such as significantly higher or lower coverage) may indicate folding or expanded repeats respectively. (2) The presence of a large number of insertions and deletions (indels) in certain regions may indicate assembly errors because correctly assembled contigs should show consistency with reads during alignment. (3) Paired reads from the same double-stranded DNA molecule should be aligned within a reasonable distance and have the correct direction on the same contig. Incorrect assembly may lead to abnormal read alignments, such as paired reads being aligned to different contigs, unreasonable insertion sizes of paired reads aligned to the same contig, or incorrect paired directions of reads. Therefore, the present invention first extracts multiple characteristic signals of the candidate error region.
[0023] 1. Regarding coverage anomalies Generally, reads can be uniquely aligned to a contig, while reads of repetitive sequences may have multiple alignment positions. Coverage is the number of reads aligned at a position. Correct assembly is usually accompanied by uniform coverage with a very small fluctuation range. If there are folded or expanded repeats in the candidate error region, abnormal coverage higher or lower than the average coverage of adjacent regions may occur.
[0024] For this reason, extract the coverage anomaly signal of the candidate error region, specifically: Using the candidate error region ( ) as a reference, extend a second preset length (such as 400bp) to both sides of it to form an extended region ( , and the total length of the extended region is 1000bp); Calculate the average coverage of all sites within the extended region ( ), traverse and detect the coverage of each site within the candidate error region. If the coverage ( ) of a certain site i is not less than the first preset multiple (such as 2 times) of the average coverage, or not greater than the second preset multiple (such as 0.5 times) of the average coverage, then this site is marked as a coverage anomaly site (i.e., site type Marked as 1); Statistically count the total number of sites with abnormal coverage within the candidate error region, calculate the ratio of the total number to the total length of the candidate error region, and use this ratio as the coverage abnormal signal of the candidate error region ( ).
[0025] It can be achieved through the following formula: ; ; ; ; ; In the formula, , , , are the start position of the candidate error region, the end position of the candidate error region, the start position of the extended region, and the end position of the extended region respectively; i refers to a specific base site in the region to be detected. According to the BAM file, the alignment information in the region to be detected is preprocessed and saved in the array Base. Then, the corresponding information is extracted from the Base array base by base in the region to be detected and the corresponding abnormal signal is calculated.
[0026] 2. Regarding insertions or deletions (Indels) anomalies There may be a large number of insertions or deletions in the candidate error region. Insertion and deletion events mean that the aligned reads show insertions or deletions relative to the contig. For example, if the contig has the sequence CGCA, but the aligned reads at the corresponding position indicate an inserted TG sequence (CTGGCA), then the contig may contain a deletion error.
[0027] Therefore, extract the insertion or deletion abnormal signal of the candidate error region, specifically: Statistically count the number of reads with insertions or deletions at each site in the candidate error region. If the coverage of a certain site is not less than the preset value, then calculate the ratio of the number of reads ( ), to the coverage of this site ( ), and use this ratio as the insertion or deletion signal of this site ( ); Compare the insertion or deletion signals of all sites in the candidate error region, and use the maximum ratio as the insertion or deletion abnormal signal of the candidate error region (
[0028] It can be achieved through the following formula: ; ; Among them, for a certain low-coverage site, when its coverage is less than a preset value, such as less than 0.2 times the average coverage, it indicates that the number of effective reads at this position is insufficient, which may lead to insufficient analysis results and further increase the risk of false positives. To avoid misjudgment, such sites are not considered in subsequent analyses.
[0029] In the actual assembly situation, the interference of heterozygous genes in diploid organisms may also generate certain indel signals, resulting in false positives. At a certain site on a contig, when there are a large number of insertion or deletion events in the corresponding reads, these reads can be identified as "error-supporting" reads (Indel-supporting Reads), that is, reads that show significant indel differences from the contig in the alignment. Compared with genetic variations, the expected proportion of "error-supporting" reads in assembly errors is higher ( ratio indel > 0.5). By setting a threshold, filtering operations are performed on the insertion-deletion signals ( ), and regions with low indel intensity are filtered out and do not participate in subsequent clustering analyses to reduce the occurrence of false positives.
[0030] 3. Abnormal pairing direction When reads are aligned with a contig, the reads should be paired in the correct direction. Mismatching of different strands will result in a large number of reads with incorrect directions in the alignment, and these reads usually gather in a local area of the candidate error region.
[0031] Therefore, the abnormal pairing direction signals in the candidate error region are extracted, specifically: Detect the reads with mismatched directions at each site in the candidate error region one by one. If reads with mismatched directions are found at a certain site, the local area of this site is marked as the error pairing direction region, and the number of reads with mismatched directions in the error pairing direction region and the total number of reads in the error pairing direction region are counted. The ratio of the number of reads ( ) to the total number ( ) is used as the abnormal pairing direction signal ( ) of this error pairing direction region; Compare the sizes of the abnormal pairing direction signals of all error pairing direction regions within the candidate error region. Preferably, the "error pairing direction regions" where the proportion of all reads with abnormal pairing directions is greater than 20% are selected, and the abnormal pairing direction signal of the error pairing direction region with the maximum value is used as the abnormal pairing direction signal ( ) of this candidate error region.
[0032] It can be achieved through the following formula: ; .
[0033] 4. Regarding abnormal insert size The insert size refers to the distance between the two paired ends of paired reads. Under normal circumstances, this distance should fluctuate within a certain range. If the insert size is too short or too long, the insert size is unreasonable. The distance of the insert size of reads can be considered to follow a normal distribution. If the insert size of reads is abnormal and deviates from the limit of normal reads (i.e., less than the minimum value of the insert size of normal reads or greater than the maximum value of the insert size of normal reads ), assembly errors may occur.
[0034] To this end, extract the abnormal insert size signal of the candidate error region, specifically: Detect the reads with abnormal insert size at each site in the candidate error region one by one. If reads with abnormal insert size are found at a certain site, mark the local area of this site as the abnormal insert size region, count the number of reads with abnormal insert size in the abnormal insert size region and the total number of reads in the abnormal insert size region, and use the ratio of the number of reads ( ) to the total number ( ) as the abnormal insert size signal ( ) of this abnormal insert size region; Compare the sizes of the abnormal insert size signals of all abnormal insert size regions within the candidate error region. Preferably, select the "abnormal insert size region" where the proportion of all reads with abnormal insert size is greater than 20%, and use the abnormal insert size signal of the abnormal insert size region with the maximum value as the abnormal insert size signal ( ) of this candidate error region.
[0035] It can be achieved through the following formula: ; ; ; ; In the formula, , are the mean and standard deviation of this sequencing data respectively.
[0036] 5. Regarding abnormal read pairing Paired reads of the same DNA molecule usually align to the same contig, and incorrect assembly of the contig may lead to a large number of orphan reads aligning to different contigs.
[0037] Therefore, extract the read pair abnormal signals in the candidate error regions, specifically: Detect the orphan reads in the candidate error regions site by site. If paired reads that cannot be aligned to the same assembly result are found at a certain site, mark the local region of this site as the orphan read region, count the number of abnormally paired reads in the orphan read region and the total number of reads in the orphan read region, and use the ratio of the number of reads ( ) to the total number ( ) as the read pair abnormal signal in the orphan read region ( ); Compare the magnitudes of the read pair abnormal signals in all orphan read regions within the candidate error region. Preferably, select the "orphan read regions" where the proportion of all abnormally paired reads is greater than 20%, and use the read pair abnormal signal with the maximum value in the orphan read region as the read pair abnormal signal in this candidate error region ( ).
[0038] It can be achieved through the following formula: ; .
[0039] 6. Regarding clipping anomalies During the alignment process, reads of some sequences that cannot be fully aligned to the contig usually occur in complex regions of the genome, such as insertion, deletion, or duplication regions. Clipping means that when the aligned reads are being aligned, one or both ends of them have unaligned sequences.
[0040] Therefore, extract the clipping abnormal signals in the candidate error regions, specifically: Count the number of reads with clipping in the candidate error regions site by site. If the coverage at a certain site is not less than the preset value, calculate the ratio of the number of reads ( ) to the coverage at this site ( ), and use this ratio as the clipping abnormal signal at this site ( ); Compare the clipping abnormal signals of all sites in the candidate error region, and use the maximum ratio as the clipping abnormal signal in this candidate error region ( ).
[0041] It can be achieved through the following formula: ; ; Among them, for a low-coverage site, when its coverage is less than a preset value, such as less than 0.2 times the average coverage, it may lead to untrustworthy analysis results, thereby increasing the false positive risk. To avoid misjudgment, such sites are not considered in subsequent analyses.
[0042] S2: For each candidate error region, based on multiple characteristic signals extracted from the candidate error region, after clustering, the clustering results of each characteristic signal of the corresponding candidate error region are obtained; if the clustering result of one characteristic signal is abnormal, then this candidate error region is an incorrect assembly region.
[0043] In the case where all candidate error regions are assembly errors, since clustering such as K-means may misclassify some real incorrect assemblies as normal. To avoid this problem, before clustering, it also includes: Randomly select several regions from the assembly results longer than the preset sequence length and jointly perform abnormal signal clustering analysis with the candidate error regions. A certain number of normal regions are introduced into the candidate dataset to assist clustering, so as to enhance the discrimination ability of incorrect assembly detection.
[0044] Clustering can group objects with similar characteristics together. K-means clustering aims to divide into two categories: normal regions and incorrect regions according to the abnormal signals extracted from the candidate error regions. K-means performs well in low dimensions, and the selected abnormal signals are independent of each other. Therefore, each candidate error region is separately clustered according to its multiple signal characteristics (such as abnormal coverage, indels, incorrect pairing direction, etc.) and the clustering results of the corresponding characteristic signals are output.
[0045] Analyze the clustering results of multiple abnormal signal characteristics of the candidate error regions respectively, and the abnormal characteristic information extracted from each region, and conduct a comprehensive evaluation of the assembly errors of multiple signals. For each candidate error region of the genome, independent clustering results of six abnormal signals (coverage, indels, pairing direction, insert size, paired reads, and clipping) of this region are obtained { , , , , , , where 0 represents normal and 1 represents abnormal. If the clustering results of these abnormal signals are all normal, it is considered that there is no abnormality in this region; compared with the signals in the correctly assembled region, the signals in the misassembled region show significant characteristic differences. Clusters containing real assembly errors often show significant performance in a certain signal feature. If at least one abnormal signal is clustered as abnormal, it is determined that there is an assembly error in this region. By defining an overall evaluation function Fregion to summarize the K-means clustering results of multiple signals and give the final abnormality evaluation, the accuracy of assembly error detection is improved, and the missed detection or misjudgment caused by the misjudgment of a single signal is avoided.
[0046] 。
[0047] S3: Based on the misassembled region, determine the type of assembly error. If the type is insertion or deletion, remove the inserted sequence or the inserted and deleted fragment; If the type is misconnection, generate a consensus sequence through multiple sequence alignment, compare the consensus sequence with the assembly result, and perform repair.
[0048] Problems such as insertions, deletions, or misjoins may exist in the assembled genome. These assemblies will disrupt the gene structure and thus affect downstream functional analysis. Therefore, after accurately locating the misassembled region, different repair strategies are adopted for different error types to improve the consistency and accuracy of the assembled sequence.
[0049] Insertion or deletion errors usually show significant differences from the reads alignment signals. In the BAM file, extract the reads sequences supporting the insertion or deletion by parsing the CIGAR signal, calculate the proportion of reads supporting this error, and mark the abnormal region. Statistically calculate the support ratio of the reads supporting the assembly error. If the support ratio at a certain position exceeds the preset threshold (such as >0.5), then perform the correction operation, that is, for insertion errors, remove the redundant sequence in the contig; for deletion errors, insert the missing fragment to restore integrity. If multiple indel error signals are detected and their support frequencies do not exceed the threshold, no indel error repair is performed.
[0050] Misjoin usually occurs during the assembly process. Common forms include cross-chromosome splicing, reverse splicing, or overly long-distance splicing. These errors may be caused by the assembly algorithm or sequence complexity, and are manifested as abnormal contig breakpoints, alignment directions, or overly long splicing distances in the alignment results.
[0051] To accurately identify the misconnected regions, first analyze and compare the characteristic signals (such as an increase in clipping events, a drastic change in coverage, etc.), and mark the breakpoint regions (1.5 kb upstream and downstream are used as the analysis range); Subsequently, based on the paired-end sequencing data, extract the paired-end reads that completely cover the breakpoint regions. After filtering, obtain the high-precision local read sequences; Use the SPOA tool to perform multiple sequence alignment on the local read sequences, generate a consensus sequence, and evaluate its integrity and consistency. Compare the consensus sequence with the assembly result. If the successfully aligned length of the repaired sequence accounts for more than 90% of the breakpoint region, it is considered that the repair is effective. Use the generated PAF file from the alignment to replace the misconnected region in the assembly result with the consensus sequence to ensure the continuity and accuracy of the overall sequence.
Claims
1. A method for detecting and repairing error assembly based on double-end sequencing data, characterized in that: The following steps are involved: S1: Based on the assembly results, the preset suspicious site is taken as the core and extended to both sides of the genome sequence by the first preset length as the candidate error region; The double-end sequencing data is compared with the assembly results, and based on the comparison information, multiple characteristic signals of the candidate error region are extracted, wherein the multiple characteristic signals are coverage abnormality signal, insertion or deletion abnormality signal, pairing direction abnormality signal, insertion size abnormality signal, read pairing abnormality signal, and shearing abnormality signal; S2: For each candidate error region, clustering is performed based on the multiple characteristic signals extracted from the candidate error region to obtain the clustering result of each characteristic signal of the corresponding candidate error region; if the clustering result of one characteristic signal is abnormal, the candidate error region is an incorrect assembly region; S3: Based on the incorrectly assembled region, determine the type of assembly error. If the type is insertion or deletion, remove the inserted sequence or the indel fragment; If the type is incorrect connection, a consensus sequence is generated through multiple sequence alignment, and the consensus sequence is aligned with the assembly result for repair.
2. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: S1, extracting the coverage anomaly signal of the candidate error area, specifically: Taking the candidate error area as a reference, extending the candidate error area to both sides of the candidate error area by a second preset length to form an extended area; Calculate the average coverage of all sites in the extended area, traverse the coverage of each site in the candidate error area, and if the coverage of a site is not less than the first preset multiple of the average coverage, or not greater than the second preset multiple of the average coverage, then the site is recorded as a coverage abnormal site; The total number of coverage anomaly sites in the candidate error region is counted, and the ratio of the total number to the total length of the candidate error region is calculated, and the ratio is used as the coverage anomaly signal.
3. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: S1, extracting insertion or deletion abnormal signals of candidate error regions, specifically: The number of reads with insertion or deletion in the candidate error region is counted site by site. If the coverage of a site is not less than the preset value, the ratio of the number of reads to the coverage of the site is calculated and the ratio is used as the insertion or deletion signal of the site. The insertion or deletion signals of all sites in the candidate error region are compared, and the maximum ratio is taken as the insertion or deletion abnormal signal of the candidate error region.
4. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: S1, extracting the pairing direction abnormal signal of the candidate error area, is specifically: Detect the reads with mismatched directions in the candidate error regions one by one. If a read with mismatched directions is found in a certain site, the local area of the site is recorded as the wrong pairing direction region. The number of reads with mismatched directions in the wrong pairing direction region and the total number of reads in the wrong pairing direction region are counted, and the ratio of the number of reads to the total number is taken as the pairing direction abnormality signal of the wrong pairing direction region. The sizes of the pairing direction abnormal signals of all the wrong pairing direction regions in the candidate error region are compared, and the pairing direction abnormal signal of the wrong pairing direction region with the maximum value is taken as the pairing direction abnormal signal of the candidate error region.
5. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: The S1, extracting the insertion size abnormal signal of the candidate error area, is specifically: Detect reads with abnormal insertion size in candidate error regions site by site. If reads with abnormal insertion size are found in a certain site, the local area of the site is recorded as the abnormal insertion size region. The number of reads with abnormal insertion size in the abnormal insertion size region and the total number of reads in the abnormal insertion size region are counted, and the ratio of the number of reads to the total number is used as the abnormal insertion size signal of the abnormal insertion size region. The insertion size anomaly signal sizes of all the insertion size anomaly regions in the candidate error region are compared, and the insertion size anomaly signal of the insertion size anomaly region with the maximum value is used as the insertion size anomaly signal of the candidate error region.
6. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: S1, extracting read segment pairing abnormality signals of candidate error regions, specifically: Detect isolated reads in candidate error regions site by site. If a site is found to be unable to match paired reads of the same assembly result, the local area of the site is recorded as an isolated read region. The number of reads with paired abnormalities in the isolated read region and the total number of reads in the isolated read region are counted, and the ratio of the number of reads to the total number is used as the signal of abnormal read pairing in the isolated read region. The read pairing anomaly signals of all isolated read regions in the candidate error region are compared, and the read pairing anomaly signal of the isolated read region with the maximum value is taken as the read pairing anomaly signal of the candidate error region.
7. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: The S1, extracting the shearing abnormality signal of the candidate error area, is specifically: The number of reads that have been sheared in the candidate error region is counted site by site. If the coverage of a site is not less than the preset value, the ratio of the number of reads to the coverage of the site is calculated and the ratio is used as the shearing abnormality signal of the site. The shearing abnormality signals of all sites in the candidate error region are compared, and the maximum ratio is taken as the shearing abnormality signal of the candidate error region.
8. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: Said S2, before clustering, also includes: Several regions are randomly selected from the assembly results with a sequence length greater than a preset length, and abnormal signal clustering analysis is performed together with the candidate error regions.
9. The method for detecting and repairing error assembly based on double-end sequencing data according to claim 1, characterized in that: In S3, if the type is wrong connection, a consensus sequence is generated through multiple sequence alignment, and the consensus sequence is aligned with the assembly result for repair, specifically: Based on the assembly results of the wrong connection type, analyze and align the characteristic signals and mark the breakpoint regions; Based on the paired-end sequencing data, the paired-end reads that completely cover the breakpoint region are extracted and filtered to obtain the local read sequence; The local read sequence is aligned through multiple sequences to generate a consensus sequence, which is then aligned with the assembly result. If the repair sequence successfully aligns for more than 90% of the breakpoint region, the consensus sequence is used to replace the incorrectly connected region in the assembly result.
10. A system for detecting and repairing error assembly based on double-end sequencing data, characterized in that: include: Feature signal extraction module: Based on the assembly results, the module takes the preset suspicious sites as the core and extends the first preset length to both sides of the genome sequence as candidate error regions; The double-end sequencing data is compared with the assembly results, and based on the comparison information, multiple characteristic signals of the candidate error region are extracted, wherein the multiple characteristic signals are coverage abnormality signal, insertion or deletion abnormality signal, pairing direction abnormality signal, insertion size abnormality signal, read pairing abnormality signal, and shearing abnormality signal; Error assembly region determination module: for each candidate error region, clustering is performed based on the multiple characteristic signals extracted from the candidate error region to obtain the clustering result of each characteristic signal of the corresponding candidate error region; if the clustering result of one characteristic signal is abnormal, the candidate error region is the error assembly region; Repair module: Based on the incorrectly assembled region, determine the type of assembly error. If the type is insertion or deletion, remove the inserted sequence or the inserted and deleted fragments; If the type is incorrect connection, a consensus sequence is generated through multiple sequence alignment, and the consensus sequence is aligned with the assembly result for repair.
Citation Information
Patent Citations
Method and device for repairing genome sequencing and assembling results, and storage medium
CN110310702A
Detection and correction system based on metagenome splicing error
CN114155914A
Whole genome structure variation identification method based on three-generation sequencing
CN115831222A
Method and system for identifying genomic sequence classification errors based on machine learning
CN115910216A
Method and system for detecting shear interval variation based on allele perception
CN119495356A