A computing device and method for identifying artificial chimeric reads in FFPE sample second-generation sequencing data

Abstract

The present invention discloses a computing device and method for identifying artificial chimeric reads in the second generation sequencing data of FFPE samples. The computing device includes three modules: a construction module of an artificial chimeric read database, an identification module of false positive mutations in a VCF file, and a correction module of artificial chimeric reads in an alignment file. Using a computing device to identify artificial chimeric reads in FFPE samples includes two modes: a fast calculation mode and a comprehensive correction mode. The fast calculation mode includes a construction module of an artificial chimeric read database and an identification module of false positive mutations in a VCF file. The comprehensive correction mode includes a construction module of an artificial chimeric read database and a correction module of artificial chimeric reads in an alignment file. Different modes are selected for different gene detection ranges. The present invention can effectively reduce the false positive rate, eliminate false positive mutations in FFPE samples caused by enzyme cutting, and reduce interference with subsequent clinical research, while ensuring that the true positive rate is not affected.

Description

Technical Field

[0001] The present invention belongs to the field of gene detection technology, and in particular relates to a computing device and method for identifying artificial chimeric reads in second-generation sequencing data of FFPE samples. Background Art

[0002] FFPE (formalin-fixed paraffin-embedded) samples are the most commonly used tissue sample preparation method in pathology laboratories. Nucleic acid molecules from FFPE are often used as retrospective and prospective research materials for clinical research. With the rapid development of NGS (next-generation sequencing), sequencing costs have decreased exponentially, and more and more NGS data from FFPE samples are used for clinical analysis. However, after treating tissue samples with formalin and paraffin, nucleic acid molecules will be significantly degraded, damaged, or chemically modified. Damage to nucleic acid molecules can be treated by DNA (Deoxyribonucleic acid) repair enzymes, but single-stranded molecules formed by DNA degradation may cause false positive mutations in downstream analysis.

[0003] DNA library preparation is a key step in NGS, and one of the more important steps is to shear DNA to produce DNA fragments of uniform length. Ultrasonic shearing is a commonly used shearing method, but in addition to the one-time large investment in purchasing an ultrasonic shearing instrument, the shearing tubes required for each shearing are also a relatively large continuous investment. In addition, the ultrasonic shearing method can only process a maximum of a dozen samples at the same time with the same instrument. For laboratories with higher throughput, it is necessary to invest in the purchase cost of multiple shearing instruments and high labor costs. Therefore, many research groups and reagent companies have developed an alternative to ultrasonic shearing: enzyme shearing. Enzyme shearing is a method of mixed enzymes, using restriction endonucleases to recognize and cut the chemical bond between phosphate and deoxyribose. For DNA enzyme shearing, it is necessary to control the amount of enzyme input and the duration of enzyme shearing. Restriction endonucleases are usually located in bacteria and are used to defend against the invasion of exogenous DNA. Restriction endonucleases can cut DNA into fragments of different lengths, producing smooth ends or sticky ends. DNA enzyme digestion is easy to use, low cost, highly scalable, and can be used in 96-well plates to process large quantities of samples at one time. Therefore, in NGS molecular laboratories, enzyme digestion can be used for automated library construction, which can process large quantities of samples at one time, thereby reducing labor costs and improving experimental efficiency.

[0004] However, DNA enzyme digestion has systemic bias and can cause primer sequencing errors. In the treatment of FFPE samples using enzyme digestion, many unexpected SNVs (Single-nucleotide variants) and INDELs (insertion or deletion), i.e., false positive mutations, were found. These mutations come from single-stranded DNA sequences in different regions of FFPE degradation, and are combined through partial matching of the end sequences to form artificial chimeric sequences, which are specifically manifested in the IGV (Integrative Genomics Viewer) image showing that there are more base sequences inconsistent with the reference gene at the 5' or 3' end of the read. This type of false positive mutation often has more read support, and the mutation frequency is generally below 15%. When studying tumor somatic mutations in FFPE samples, since the lower limit of mutation frequency detection of tumor somatic mutations is generally required to be above 5%, these false positive mutations will affect subsequent clinical research.

[0005] Currently, there is a method that simply filters by removing mutations within 15bp (base pair) from both ends of the read. However, due to the specificity of FFPE samples and inconsistent mixed enzymes, this method still retains a large number of false positive mutations and may cause the loss of true positive mutations. Summary of the invention

[0006] In view of the situation that FFPE sample enzyme cutting and interruption have artificial chimeric reads leading to many false positive mutations, the present invention provides a computing device and method for identifying artificial chimeric reads in the second-generation sequencing data of FFPE samples, which can effectively reduce the false positive rate and eliminate the false positive mutations of FFPE samples caused by enzyme cutting while ensuring that the true positive rate is not affected, reduce the interference with subsequent clinical research, and can reflect the quality of FFPE samples through three indicators of IPercent, FPercent or APercent, which are used to evaluate the reliability of the experiment.

[0007] In order to achieve the above purpose, this application adopts the following technical solutions:

[0008] In a first aspect, the present invention provides a computing device for identifying artificial chimeric reads in next-generation sequencing data of FFPE samples, the computing device comprising three modules, namely, a module for constructing an artificial chimeric read database, a module for identifying false positive mutations in a VCF file, and a module for correcting artificial chimeric reads in an alignment file;

[0009] The computational device can be used to identify artificial chimeric reads in FFPE samples in two modes: fast computational mode and comprehensive correction mode.

[0010] The fast calculation mode includes a module for building an artificial chimeric read database and a module for identifying false positive mutations in VCF files.

[0011] The comprehensive correction mode includes a construction module for the artificial chimeric read database and a correction module for artificial chimeric reads in the alignment file;

[0012] Different modes are selected for different gene detection ranges. For small probe combinations, the comprehensive correction mode is selected, and for whole exome sequencing or large probe combinations, the fast calculation mode is selected.

[0013] In the above technical solution, the construction module of the artificial chimeric read database is divided into two parallel computing parts: artificial chimeric sequence cluster region identification and artificial chimeric read database construction, including the following contents:

[0014] S11, sequencing data is generated from FFPE samples after nucleic acid extraction, enzyme shearing, library preparation, hybridization capture and sequencing;

[0015] S12, the sequencing data is aligned with the hg19 reference genome by bwa to obtain an alignment file, and the repeated fragments in the alignment file are marked by picard to obtain an alignment file with marked repeated fragments;

[0016] S13, divide the potential parallel computing units according to different computing modes:

[0017] S131, if it is a fast computing mode, divide the potential parallel computing units according to the mutations retained in the VCF, and the specific operations are as follows:

[0018] (1) Retain relatively realistic candidate detection mutations, and retain mutations that meet the following conditions: the filter information FILTER is PASS, the mutation fragment support number AD is greater than or equal to a certain value, the site depth DP is greater than or equal to a certain value, and the mutation frequency AF is greater than or equal to a certain value;

[0019] (2) The region where the above-mentioned retained mutation is located is used as a potential parallel computing unit. Assuming that the chromosome number of the mutation is chr, the mutation position is pos, the reference base is ref, and the mutation base is alt, then the chromosome number of the potential parallel computing unit is chr, the starting position number is pos-1, and the ending position number is pos-1 plus the length of ref. The mutation information is stored in the structure of chr-pos-ref-alt;

[0020] S132, if it is a comprehensive correction mode, capture the sequencing probe combination sequencing for the target region, and use all the target regions as potential parallel computing units;

[0021] S14, merging the above potential parallel computing units according to the following method to reasonably reduce the overall running time:

[0022] (1) Assume that a potential parallel computing unit i, chromosome is chr i , the starting position is start i , the end position is end i , another potential parallel computing unit j, chromosome is chr j , the starting position is start j , the end position is end j , if the above two units meet any of the following conditions, the two calculation units are merged:

[0023] a、chr i =chr j And start i -end j ≤PE;

[0024] b、chr i =chr j And start j -end i ≤PE;

[0025] The above PE is the read length of paired-end sequencing;

[0026] (2) If the number of potential computing units Num corresponding to the merged computing unit exceeds a certain value, the computing unit is split according to the coordinate position sorting to ensure that Num has a certain upper limit and the running time of each merged computing unit is relatively balanced. If it is a fast computing mode, for the parallel computing unit, the mutation information corresponding to the unit is stored at the same time;

[0027] S15, parallel processing for each computational unit to find all artificial chimeric sequence clusters in the region;

[0028] S16, the artificial chimeric sequence clusters are numbered according to chromosome number chr i And the genome alignment position pos i Merge into a region set AR containing multiple artificial chimeric sequence clusters = {AR k ,k∈T *},

[0029] Where k represents the region after the kth artificial sequence cluster is merged, T * Represents the total number of all artificial chimeric sequence regions in the FFPE sample, and each artificial chimeric sequence cluster region AR k ={AC i ,i∈R k}, AC i Represents an artificial chimeric cluster sequence in the region, stores sequence-related information, i represents the i-th artificial chimeric sequence cluster, R k The number of artificial chimeric sequence clusters contained in the corresponding k-th region,

[0030] Total number of artificial chimeric cluster sequences

[0031] S17, program design is performed according to the monitoring mode, a computing unit is specifically used for monitoring queue information, and all the artificial chimeric sequence cluster regions are divided according to the number of processes to form multiple parallel artificial chimeric read segment identification computing units;

[0032] S18, for each artificial chimeric read segment identification unit, complete the following calculation process:

[0033] S181, traverse and match to chromosome number chr k , the genome alignment position is at start k and end k All reads between n and n1 were filtered out: reads that were not aligned to the reference genome, secondary aligned, supplementary aligned, marked as duplicates, or had an alignment quality less than n1 (i.e., MIN_MQ);

[0034] S182, a matrix M storing the fully matched read base positions and reference genome positions by using the read cigar value, MD label, and read alignment start position. e ;

[0035] S183, judging according to the read segment name, if the read segment is located in the artificial chimeric cluster library corresponding to the artificial chimeric region, the following judgment is performed and the corresponding information is put into the monitoring queue:

[0036] (1) If a read is in more than two artificial chimeric clusters, and the clustering direction tag ctag of the corresponding artificial chimeric cluster contains both “+” and “-”, the read is defined as a non-aligned sequence, and the read name, tag rtag (whether the read segment is read1 or read2), and the read position of the chimeric site are put into the monitoring queue;

[0037] (2) If the ctag of the read corresponding to the artificial sequence cluster is "+", then through the matrix M eGet the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the minimum value as the breakpoint position of the read, and determine the ratio of the length from the position to the 3' end of the read to the entire read length. If it is greater than or equal to a certain value, it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number chr k Put it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0038] (3) If the ctag of the read corresponding to the artificial sequence cluster is “-”, then the matrix M e Get the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the maximum value as the breakpoint position of the read. Determine the ratio of the length from the 5' end of the read to this position to the entire read length. If it is greater than or equal to a certain value, it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, tag rtag and chromosome number chr k Put it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0039] S184, for matrix M e , filter the submatrix FM containing the breakpoint information of the artificial chimeric cluster set e , if the number of submatrix rows is greater than 0, further determine whether the read is an artificial chimeric read:

[0040] (1) If the ctag of the artificial chimeric cluster is "+", check whether the length L of the sequence from the read position of the read breakpoint to the 3' end of the read is greater than or equal to n2 (i.e., CLIPLEN) and whether the sequence is equal to the cluster identification sequence CS i A partial base sequence of length L starting from the 5' end;

[0041] (2) If the ctag of the artificial chimeric cluster is “-”, check whether the length L of the sequence from the 5’ end of the read to the breakpoint is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence CS i A partial sequence of length L starting from the 3' end of

[0042] If the above conditions are met, the read is considered an artificial chimeric read;

[0043] S185, if after screening, the read meets the artificial chimeric cluster judgment conditions, the information is stored according to the three modes (1), (2), and (3) of step S183;

[0044] S19, extract information from the monitoring queue, use the chromosome number as the first-level index, and the read name and read rtag as the second-level index to build an artificial chimeric read database; for different data, store the data corresponding to the index according to the following method:

[0045] For the queue information whose read name is consistent with the read rtag, if there is information marked as non-matching, the read is added with the non-matching tag and stored in the database;

[0046] If the ctag is both "+" and "-", the read will also be tagged with a non-matching tag and stored in the database;

[0047] If only "+" or "-" exists, the queue information is sorted from small to large according to the read position of the chimeric site. If the ctag is "+", the first information is stored. If the ctag is "-", the label ctag of the last queue information and the read position of the chimeric site are stored.

[0048] In the above technical solution, the recognition module for false positive mutations in the VCF file constructs a process pool according to the artificial sequence cluster recognition units divided in the construction module of the artificial chimeric read database, and uses parallel computing to determine whether a mutation is a false positive mutation caused by an artificial chimeric sequence;

[0049] If it is caused by artificial chimeric sequences, the VCF information will be modified and the new information will be stored in the result file. If there is no impact, the original information will be retained in the result file;

[0050] The specific steps include:

[0051] S21, performing parallel calculations according to the calculation units divided in the construction module of the artificial chimeric read database to obtain a P value of whether each mutation is affected by the chimeric sequence and various data of mutations excluding the influence of the chimeric sequence:

[0052] (1) Assume that the unit is calculated as unit i, and the chromosome is chr i , the starting position is start i , the end position is end i , including M i mutation sites;

[0053] (2) Chromosome number chr iAs an index, query the artificial chimeric read database to store the chimeric read information set corresponding to the chromosome;

[0054] (3) Construct an M i ×10 matrix Matrx i Used to store various types of data for each mutation site, including: the number of reads that support the ref type and are read1, the number of reads that support the ref type and are read2, the number of reads that support the ref type and are read1 and are not artificial chimeric reads, the number of reads that support the ref type and are read2 and are not artificial chimeric reads, the number of reads that support the alt type and are read1, the number of reads that support the alt type and are read2, the number of reads that support the alt type and are read1 and are not artificial chimeric reads, the number of reads that support the alt type and are read2 and are not artificial chimeric reads, the sequencing depth of the site, and the depth of the site excluding artificial chimeric sequences;

[0055] (4) Iterate the chromosome number. The chromosome is chr i , starting position start i and end position end i All reads between the two sequences were filtered out: those that were not aligned to the reference genome, those that were secondary aligned, those that were supplementally aligned, those that were marked as duplicates, or those that had an alignment quality less than MIN_MQ;

[0056] (5) For the retained reads, use the matrix Matrx i Complete multi-threaded judgment, first for each row mutation m (m∈M i ), the mutation information is chr-pos-ref-alt, assuming that the first base of ref and alt is the same, the starting position of the mutation altstart is set to pos-1, if not, altstart is set to pos-2, and the end position of the mutation altend is set to pos plus the length of ref;

[0057] (6) Whether the partial sequence of the read at the mutation position is a chimeric sequence, the read positions corresponding to altstart and altend are obtained by using the cigar value and the start position of the read chromosome alignment: readstart and readend; if the read name plus rtag exists in the artificial chimeric read library and meets the following conditions, it is marked as a chimeric sequence:

[0058] a. Reads are marked as non-aligned;

[0059] b. The ctag of the Read is "+", and the read position rpos of the chimeric breakpoint is less than or equal to readstart;

[0060] c. The ctag of the Read is "-", and the read position rpos of the chimeric breakpoint is greater than or equal to readend;

[0061] (7)Matrx i The mth row and the 9th column are incremented by 1. If it is judged not to be a chimeric sequence, the mth row and the 10th column are incremented by 1;

[0062] (8) Determine whether the sequence between readstart and readend of read is consistent with ref or alt. If it is consistent with ref, it is marked as isref; if it is consistent with alt, it is marked as isalt.

[0063] (9) If read is read1, set the column position col to 1. If read is read2, set col to 2. If read is marked as isalt, add 4 to the col position. The matrix Matrx i Add 1 to the mth row and colth column;

[0064] (10) If it is not a chimeric sequence, add 2 to col and the matrix Matrx i Add 1 to the mth row and colth column;

[0065] (11) Through the above method, all reads in the calculation unit i area are traversed to obtain the matrix Matrx i ;

[0066] (12) For each mutation, assume that the mutation is located in Matrx i If there are m rows in the matrix, a new two-dimensional matrix Data is constructed, and its calculation rules are as follows:

[0067]

[0068] (13) Data was used as input data and one-sided Fisher's exact test was used to calculate the oddratio and P value;

[0069] (14) If oddratio is less than a certain value, the mutation is considered to be affected by the artificial chimeric sequence. Under this condition, if the P value is less than a certain value, the mutation is considered to be significantly affected by the artificial chimera.

[0070] (15) Construct storage information EAVar for the affected mutations i The storage information is indexed by chr-pos-ref-alt, and the value contains Matrx i [m][3], Matrx i [m][4]、Matrxi [m][7], Matrx i [m][8]、Matrx i [m]

[10] , oddratio, and P value;

[0071] S22, modify the VCF file, add the ISAE tag to the INFO part of the header file, indicate that the mutation is affected by the chimeric sequence, add the isAEM tag to the FILTER part, indicate that the mutation is a false point caused by the chimeric sequence and is filtered, and then modify the mutation information of each computing unit. Assuming there are K computing units and each computing unit runs in parallel, the modification steps are as follows:

[0072] (1) Retain relatively realistic candidate detection mutations, and retain mutations that meet the following conditions: FILTER is PASS, the AD of the mutation is greater than or equal to a certain value, the DP is greater than or equal to a certain value, and the AF is greater than or equal to a certain value;

[0073] (2) Store the total number of mutations Num retained by each computing unit i ;

[0074] (3) Determine whether the retained mutation is located in EAVar i In the database, if located, the ISAE label is added to the INFO information section to indicate that the mutation is affected by an artificial chimeric sequence;

[0075] (4) Store the total number of mutations INum affected by the artificial chimeric sequence in each computing unit i ;

[0076] (5) According to EAVar i Database information, modify the mutated AD, DP and AF information, and the calculation rules are as follows;

[0077] AD=Matrx i [m][7]+Matrx i [m][8]

[0078] DP=Matrx i [m]

[10]

[0079]

[0080] (6) Then determine whether the mutation meets the conditions in (1). If it does not meet the conditions, and the P value is less than 0.05, the FILTER information is PASS, then change the FILTER to the isAEM label;

[0081] (7) Store the total number of mutations FNum that are filtered out due to artificial chimeric sequences in each computing unit i ;

[0082] S23, merging the VCF files of each of the above calculation units to generate a new VCF file, summarizing the mutation information affected by all artificial chimeric sequences;

[0083] S24, calculate the mutation ratio IPercent affected by the artificial chimeric sequence of the sample and the false point ratio FPercent filtered out, and the calculation rules are as follows:

[0084]

[0085]

[0086] IPercent and FPercent are used as quality control indicators to evaluate whether mutations in the next-generation sequencing of FFPE samples are affected by artificial chimeric sequences. When both quality control indicators are greater than 50%, the quality of the FFPE sample is poor.

[0087] In the above technical solution, the correction module for artificial chimeric reads in the alignment file uses an artificial chimeric read database to modify the artificial chimeric sequence information in the alignment file, change its base information to N bases, and modify the cigar value and label information to eliminate the influence of the artificial chimeric sequence on downstream analysis, including the following steps:

[0088] S31, extract the unmapped read information into a file named unmap.bam through samtools view;

[0089] S32, dividing the artificial chimeric read database in the construction module of the artificial chimeric read database into computing units according to a type of index, i.e., chromosome number, assuming that there are N computing units, and processing the alignment files in parallel:

[0090] (1) Assume that the chromosome number corresponding to a computational unit is chr i , the output comparison file is Part i .bam;

[0091] (2) In Part i The parameter information of this processing is added to the header file information of .bam for recording;

[0092] (3) Traverse and compare to chr i All reads in Part 1, if the read has secondary alignment, supplementary alignment, and duplicate tags, these reads are written to Part 2. i .bam file and skip the following processing steps;

[0093] (4) If the alignment quality of a read is greater than or equal to MIN_MQ, the number of such reads is counted, which is ReadNum i ;

[0094] (5) If the paired read of the read is in the artificial chimeric read database and is marked as non-aligned, the read is labeled with the paired read non-aligned label and then written into Part i .bam file;

[0095] (6) If the read is in the artificial chimeric read database and is marked as non-aligned, the read is labeled as non-aligned, the alignment quality is set to 0, the cigar string is set to ".", and then Part i .bam file;

[0096] (7) If the read is in the artificial chimeric read database and is not marked as non-aligned, if ctag is a "+" tag, change the read sequence from the chimeric breakpoint to the end base to N bases, and change the string corresponding to the bases from the chimeric breakpoint to the end in the cigar value to these base numbers plus the "S" tag, and then write it to Part i .bam file;

[0097] (8) If the read is located in the artificial chimeric read database and is not marked as non-aligned, if ctag is the "-" label, the bases from the beginning of the read to the chimeric breakpoint are changed to N bases, and the reference genome alignment start position of the read is set to the reference genome position corresponding to the next base of the chimeric breakpoint. At the same time, the cigar value from the beginning to the chimeric breakpoint base corresponding to the string is changed to the number of these bases plus the "S" label, and then written to Part i .bam file;

[0098] (9) The number of reads in (6), (7), and (8) is AReadNum i ;

[0099] S33. Calculate the percentage of chimeric reads in the sample, APercent, as an indicator for sample quality control evaluation. The calculation method is as follows:

[0100]

[0101] S34. Merge unmap.bam file and all Parts i .bam file to the same output comparison file out.bam, and set the header file of the output file to Part iThe header file of the .bam file, use samatools sort to sort the output comparison file out.bam, and use samtools index to build the index file of the output file;

[0102] S35. The sorted output alignment file and index file are used in subsequent downstream analysis, and the file removes the influence of chimeric sequences.

[0103] In the above technical solution, the specific steps of S15 are as follows:

[0104] S151, traverse all reads in the computational unit region, and filter out reads that are not aligned to the reference genome, are secondary aligned, are supplementally aligned, are marked as duplicates, or have an alignment quality less than MIN_MQ;

[0105] S152, the following three matrices are obtained by read alignment starting position, MD label and cigar value:

[0106] (1) Matrix M that stores the fully matched read base positions and reference genome positions e ;

[0107] (2) Matrix M storing mismatched read base positions and reference genome positions s ;

[0108] (3) A matrix M storing the cigar pattern number and the number of bases corresponding to the cigar pattern c ;

[0109] S153, artificial chimeric read candidate sequences show the following 6 patterns in the Bam file:

[0110] (1) Mode 1: There are more than CLIPLEN consecutive softclip bases at the 5' end of the read segment;

[0111] (2) Mode 2: There are more than CLIPLEN consecutive softclip bases at the 3′ end of the read segment;

[0112] (3) Mode 3: There are more than CLIPLEN continuous softclip bases at the 5' end of the read segment, and there is a SNV or INDEL mutation within n3 (i.e., ALTDIS) bases downstream;

[0113] (4) Mode 4: There are more than CLIPLEN continuous softclip bases at the 3' end of the read segment, and there is a SNV or INDEL mutation within ALTDIS bases upstream;

[0114] (5) Mode 5: There are n4 or more SNV or INDEL mutations (i.e., MINSNV, MINSNV is greater than or equal to 2) within the 5' end of the read segment ALTDIS.

[0115] (6) Mode 6: There are more than MINSNV SNV or INDEL mutations within the ALTDIS bases at the 3' end of the read segment;

[0116] S154, for the candidate chimeric reads of the six modes in S153, the following information is recorded through the three matrices in S152:

[0117] (1) Chromosome number of the reference genome to which the read segment is aligned;

[0118] (2) Obtain the read position rpos and the genome alignment position pos of the read chimeric breakpoint. The corresponding processing steps are as follows:

[0119] a. If it is mode 1, the read position of the base 1 bp downstream of the softclip sequence and the corresponding genome alignment position are recorded;

[0120] b. If it is mode 2, the read position of the base 1 bp upstream of the softclip sequence and the corresponding genomic alignment position are recorded;

[0121] c. If it is mode 3, obtain the farthest mutant base B within ALTDIS bases downstream of the 5' end softclip sequence of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0122] d. If it is mode 4, obtain the farthest mutant base B within ALTDIS bases upstream of the 3' softclip sequence of the read segment i , search upstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0123] e. If it is mode 5, obtain the farthest mutant base B within ALTDIS bases downstream of the 5' end of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j If the first base of the read segment is at B j The number of bases between the read positions is greater than CLIPLEN, then record B j The read positions and genome comparison positions of

[0124] f. If it is mode 6, obtain the farthest mutant base B within ALTDIS bases upstream of the 3' end of the read segment i , search upstream from this position to obtain the first read base B that matches the reference base j If you read segment B j The number of bases between the read position and the end of the read is greater than CLIPLEN, record B j The read positions and genome comparison positions of

[0125] (3) the name of the read segment, readname, and the label rtag indicating whether the read segment is Read1 or Read2;

[0126] (4) Clustering direction tag ctag: for modes 1, 3, and 5, the tag “-” is recorded, indicating that the subsequent clustering direction is from downstream to upstream; for modes 2, 4, and 6, the tag “+” is recorded, indicating that the subsequent clustering direction is from upstream to downstream;

[0127] (5) Chimeric sequence S: if the clustering label is “-”, the sequence from the 5’ end of the read to rpos is recorded; if it is “+”, the sequence from rpos to the 3’ end of the read is recorded;

[0128] S155, identifying candidate chimeric read clusters by clustering algorithm:

[0129] (1) Sort by chromosome chr, genome alignment position pos and cluster label;

[0130] (2) Cluster the candidate chimeric reads whose chromosome chr, genome alignment position pos and clustering tag ctag are consistent, assuming that the set R = {r i ,i∈N *}, r i Recording the candidate chimeric sequence information;

[0131] (3) The set R is sorted from large to small according to the length of the chimeric sequence S. Assume that the sorted set is RS = {rs i ,i∈N *};

[0132] (4) From rs1 to rs N Until then, iteratively calculate all candidate chimeric sequence cluster sets C, C = {C i ,i∈M *}, where C i = {rs j ,j∈K i}, The specific calculation process is as follows:

[0133] a. Assume that the candidate chimeric read of the starting iteration is rs j , the chimeric sequence is S j , the length of the chimeric sequence is L j , S j It is a set C i Identification sequence CS of candidate chimeric sequence cluster i ;

[0134] b. From rs j+1 Start, rs j+1 The chimeric sequence is S j+1 , the length of the chimeric sequence is L j+1 , if any of the following conditions are met: If the clustering direction label ctag is "-", S j+1 Equal to S j L j -L j+1 +1 to L j If ctag is "+", S j+1 Equal to S j 1st to Lth j+1 The base sequence between rs j+1 Included in set C j ;

[0135] c. If the kth rs k If the above conditions are not met, rs k The chimeric sequence S k As a new round of cluster set C i+1 The recognition sequence CS i+1 , start the next round of iteration;

[0136] (5) Finally, we get a set C, where each element C i Contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and which contains all candidate chimeric reads {rs j ,j∈K i};

[0137] S156, using local comparison to determine whether the candidate chimeric sequence cluster is a chimeric sequence, the main process is as follows:

[0138] (1) Obtain the candidate chimeric sequence cluster genome alignment position pos i Reference gene sequence RS between n5 (i.e. FLANK) bases upstream and downstream i ;

[0139] (2) The default match score is 2, the default mismatch penalty is 3, the default gap opening penalty is 10, and the default gap extension penalty is 2. The scoring matrix is ​​constructed according to the above settings;

[0140] (3) Survival cluster identification sequence CS i The reverse complementary sequence RCS i ;

[0141] (4) Sequence RCS i and reference gene sequence RS i The Smith-Waterman algorithm was used for alignment. If the sequence similarity was greater than n6 (i.e., SIM), the sequence cluster was considered to be an artificial chimeric sequence cluster of FFPE samples.

[0142] S157, according to the above method, obtain all artificial chimeric sequence cluster information sets AC = {AC i ,i∈S *}(S≤M), each artificial chimeric sequence cluster information AC i Contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and all artificial chimeric reads {rs j ,j∈K i}.

[0143] In the above technical solution, in S16, the specific merging process is as follows:

[0144] S161, sort all artificial chimeric sequence clusters from small to large according to chromosome number and genome alignment position;

[0145] S162, iterating from the first artificial chimeric cluster to the last artificial chimeric cluster, the iterative process is as follows:

[0146] (1) Assume that a certain chimeric sequence cluster is AC i , whose chromosome number is chr i , the genome alignment position is pos i ;

[0147] (2) The next human chimeric sequence cluster is AC i+1 , whose chromosome number is chr i+1 , the genome alignment position is pos i+1 , if chr i =chr i+1 And pos i +PE≥pos i+1, then AC i+1 Incorporation into AC i The artificial chimeric cluster region AR k ;

[0148] (3)AR k Contains the chromosome number chr k , genome alignment starting position start k =min(pos i ,…,pos R ), genome alignment end position endk = max(pos i ,…,pos R ) and the corresponding set of all artificial chimeric sequence clusters {AC i ,i∈R k};

[0149] (4) If an artificial chimeric sequence cluster AC j If the previous artificial chimeric sequence cluster does not meet condition B, then AC j As a new artificial chimeric region AR k+1 The next iteration starts with the element; if AR k If the number of artificial chimeric sequence clusters exceeds a certain value, the next artificial chimeric sequence cluster AC j As a new artificial chimeric region AR k+1 The next iteration starts with the element.

[0150] In a second aspect, the present invention provides a method for identifying artificial chimeric reads in the second-generation sequencing data of FFPE samples, wherein the sequencing data is generated after nucleic acid extraction, enzyme fragmentation, library preparation, hybridization capture and on-machine sequencing of the FFPE samples, and the sequencing data is processed using the above-mentioned computing device to eliminate false positive mutations caused by artificial chimeric sequences and identify artificial chimeric reads in the second-generation sequencing data of FFPE samples.

[0151] The beneficial effects of the present invention are:

[0152] 1. The present invention proposes a computing device for identifying artificial chimeric reads in the second generation sequencing data of FFPE samples, which identifies artificial chimeric sequence clusters through parallel computing units, and then constructs a database of artificial chimeric reads through parallel computing again. This database can be used in different ways to process mutation data or alignment data. By comparing with other technical platforms for detecting mutations, the TPR of the present invention can reach 100% and the FDR can reach 0%.

[0153] 2. The present invention proposes a fast computing mode for whole exome sequencing, whole genome sequencing or large probe combination sequencing technology, which reasonably divides multiple computing units and can accurately identify false positive mutations caused by artificial chimeric reads through parallel computing. The running speed for whole exome data can be as low as less than 5 minutes.

[0154] 3. The present invention proposes a comprehensive correction mode processing process, which can accurately identify chimeric sequences in the second-generation sequencing comparison data and eliminate the impact on downstream analysis by modifying the comparison data information.

[0155] 4. The present invention proposes three indicators, IPercent, FPercent and APercent, which are used to evaluate the impact of chimeric sequences from the mutation and read levels as indicators for measuring the quality of FFPE samples. BRIEF DESCRIPTION OF THE DRAWINGS

[0156] Figure 1 It is a working flow diagram of each module of the computing device of the present invention. DETAILED DESCRIPTION

[0157] In order to better illustrate the purpose, technical solutions and advantages of the present invention, the present invention will be further described below in conjunction with specific embodiments. The present invention can be implemented in many different forms and should not be construed as being limited to the embodiments set forth herein. On the contrary, these embodiments are provided so that this disclosure will be thorough and complete, and the concept of the present invention will be fully conveyed to those skilled in the art, and the present invention will only be limited by the claims.

[0158] Example 1 Fast calculation mode for whole exome sequencing data

[0159] The fast computing mode is mainly aimed at whole exome sequencing, whole genome sequencing and capturing sequencing data of target areas above the Mb (megabase, one million bases) level. This mode can complete the calculation of whole exome sequencing data within five minutes at the fastest by reasonably dividing parallel computing units and setting up multiple computing units.

[0160] (I) Construction of artificial chimeric read database

[0161] This module involves the construction of an artificial chimeric read database. This module is mainly divided into two parallel computing parts: artificial chimeric sequence cluster region identification and artificial chimeric read database construction. This module contains the following contents:

[0162] 1. The pleural and ascites samples numbered T1 and the FFPE samples numbered T2 were subjected to nucleic acid extraction, enzyme cutting, library preparation, hybridization capture and sequencing to generate sequencing data. The hybridization capture probes can be selected from probe combinations of various manufacturers, such as Aigitaikang V3 or Agilent's SureSelect Human All Exon, etc.

[0163] 2. The sequencing data was aligned with the hg19 reference genome through bwa to obtain an alignment file. The repeated fragments in the alignment file were marked by picard to obtain an alignment file with marked repeated fragments.

[0164] 3. For the alignment file of the marked duplication, use VarDict software to detect mutations in the whole exon region and obtain VCF files;

[0165] 4. Divide the potential parallel computing units. You can divide the potential parallel computing units according to the mutations retained in the VCF. The specific operations are as follows:

[0166] (1) Retain relatively realistic candidate detection mutations, and retain mutations that meet the following conditions: FILTER information is PASS, mutation AD (allele depth, mutation support number) ≥ 4 (default value, adjustable), DP (site depth, site sequencing depth) ≥ 40 (default value, adjustable), AF (allele frequency, mutation frequency) ≥ 3% (default value, adjustable);

[0167] (2) The region where the above-mentioned retained mutation is located is used as a potential parallel computing unit. Assuming that the chromosome number of the mutation is chr, the mutation position is pos, the reference base is ref, and the mutation base is alt, then the chromosome number of the potential parallel computing unit is chr, the starting position number is pos-1, and the ending position number is pos-1 plus the length of ref. The mutation information is stored in the structure of chr-pos-ref-alt;

[0168] 5. Combine the above potential parallel computing units in the following way to reasonably reduce the overall running time:

[0169] (1) Assume that a potential parallel computing unit i, chromosome is chr i , the starting position is start i , the end position is end i , another potential parallel computing unit j, chromosome is chr j , the starting position is start j , the end position is end j , if the above two units meet any of the following conditions, the two calculation units are merged:

[0170] A.chr i =chr j And start i -endj ≤PE;

[0171] B.chr i =chr j And start j -end i ≤PE;

[0172] The above PE is the read length of the paired-end sequencing. For example, if the paired-end sequencing is 150bp, the PE is 150.

[0173] (2) If the number of potential computing units (Num) corresponding to the merged computing unit exceeds a certain value, such as 50, the computing units are split according to the coordinate position to ensure that Num has a certain upper limit and the running time of each merged computing unit is relatively balanced. For the parallel computing unit, the mutation information corresponding to the unit is stored at the same time;

[0174] 6. For each computing unit, parallel processing is performed to find all artificial chimeric sequence clusters in the region. The main steps are as follows:

[0175] (1) Traverse all reads in the computational unit region and filter out reads that are not aligned to the reference genome, secondary aligned, supplementary aligned, marked as duplicates, or have an alignment quality less than MIN_MQ (the default value is 20, which can be adjusted);

[0176] (2) The following three matrices are obtained by using the read alignment starting position, MD (a set of strings that record mismatched bases and reference base missing) labels, and cigar (a description of the pattern of sequence alignment to the reference sequence and the number of bases corresponding to each type of pattern) value:

[0177] A. Matrix M storing the base positions of the reads that match perfectly (first column) and the reference genome positions (second column) e ;

[0178] B. Matrix M storing mismatched read base positions (first column) and reference genome positions (second column) s ;

[0179] C. Matrix M that stores the cigar pattern number (first column) and the number of bases corresponding to the cigar pattern (second column) c .

[0180] (3) Artificial chimeric read candidate sequences show the following six patterns in the Bam (Binary Alignment Map) file:

[0181] A. Mode 1: There are more than CLIPLEN (optional parameter, default is 6) continuous softclip bases at the 5' end of the read segment (DNA polymerase working direction, from 5' to 3' end);

[0182] B. Mode 2: There are more than CLIPLEN consecutive softclip bases at the 3' end of the read segment;

[0183] C. Mode 3: There are more than CLIPLEN continuous softclip bases at the 5' end of the read segment, and there is a SNV or INDEL mutation within the downstream ALTDIS (optional parameter, default is 30) bases;

[0184] D. Mode 4: There are more than CLIPLEN continuous softclip bases at the 3' end of the read segment, and there is a SNV or INDEL mutation within ALTDIS bases upstream;

[0185] E. Mode 5: There are more than MINSNV (optional parameter, default is 3) SNV or INDEL mutations within the 5' end of the read segment ALTDIS;

[0186] F. Mode 6: There are MINSNV or more SNV or INDEL mutations within the bases within ALTDIS at the 3' end of the read segment.

[0187] (4) For the candidate chimeric reads of the six modes in (3), the device records the following information through the three matrices in (2):

[0188] A. The reference genome chromosome number chr to which the reads are aligned;

[0189] B. Obtain the read position rpos and the genome alignment position pos of the read chimeric breakpoint. The corresponding processing steps are as follows:

[0190] a. If it is mode 1, the read position of the base 1 bp downstream of the softclip sequence and the corresponding genome alignment position are recorded;

[0191] b. If it is mode 2, the read position of the base 1 bp upstream of the softclip sequence and the corresponding genome alignment position are recorded;

[0192] c. If it is mode 3, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases downstream of the 5' softclip sequence of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0193] d. If it is mode 4, obtain the farthest mutation base (SNV or INDEL) within ALTDIS bases upstream of the 3' end softclip sequence of the read segment B i , search upstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0194] e. If it is mode 5, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases downstream of the 5' end of the read segment B i , search downstream from this position to obtain the first read base B that matches the reference base j If the first base of the read segment is at B j The number of bases between the read positions is greater than CLIPLEN, then record B j The read positions and genome comparison positions of

[0195] f. If it is mode 6, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases upstream of the 3' end of the read segment B i , search upstream from this position to obtain the first read base B that matches the reference base j If you read segment B j The number of bases between the read position and the end of the read is greater than CLIPLEN, record B j The read positions and genome comparison positions of

[0196] C. The name of the read segment readname and the tag rtag indicating whether the read segment is Read1 (sequencing sequence of the first round of template strand in double-end sequencing) or Read2 (sequencing sequence of the second round of template strand in double-end sequencing);

[0197] D. Clustering direction tag ctag, for modes 1, 3 and 5, record the tag “-”, indicating that the subsequent clustering direction is from downstream to upstream; for modes 2, 4 and 6, record the tag “+”, indicating that the subsequent clustering direction is from upstream to downstream;

[0198] E. Chimeric sequence S, if the clustering label is "-", the sequence from the 5' end of the read to rpos is recorded; if it is "+", the sequence from rpos to the 3' end of the read is recorded.

[0199] (5) Identifying candidate chimeric read clusters through clustering algorithms:

[0200] A. Prefer to sort by chromosome chr, genome alignment position pos and cluster label;

[0201] B. Cluster the candidate chimeric reads that have the same chromosome chr, genome alignment position pos, and clustering tag ctag. Suppose the set R = {r i ,i∈N *}, r i Recording the candidate chimeric sequence information;

[0202] C. The set R is sorted from large to small according to the length of the chimeric sequence S. Assume that the sorted set is RS = {rs i ,i∈N *};

[0203] D. Starting from rs1 to rs N Until then, iteratively calculate all candidate chimeric sequence cluster sets C, C = {C i ,i∈M *}, where C i = {rs j ,j∈K i}, The specific calculation process is as follows:

[0204] a. Assume that the candidate chimeric read of the starting iteration is rs j , the chimeric sequence is S j , the length of the chimeric sequence is L j , S j It is a set C i Identification sequence CS of candidate chimeric sequence cluster i ;

[0205] b. From rs j+1 Start, rs j+1 The chimeric sequence is S j+1 , the length of the chimeric sequence is L j+1 , if any of the following conditions are met: If the clustering direction label ctag is "-", S j+1 Equal to S j L j -L j+1 +1 to L j If ctag is "+", S j+1 Equal to S j 1st to Lth j+1 The base sequence between rs j+1 Included in set C j ;

[0206] c. If the kth rs k If the above conditions are not met, rs k The chimeric sequence S k As a new round of cluster set Ci+1 The recognition sequence CS i+1 , start the next round of iteration;

[0207] E. Finally, we get a set C, where each element C i Contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and which contains all candidate chimeric reads {rs j ,j∈K i};

[0208] (6) Using local comparison to determine whether the candidate chimeric sequence cluster is a chimeric sequence, the main process is as follows:

[0209] A. Obtain the candidate chimeric sequence cluster genome alignment position pos i Upstream and downstream FLANK (optional parameter, default is 500, the sum of the front and back extended sequences is 1000bp, it is generally believed that more than 1000bp may be caused by structural variation) reference gene sequence RS between bases i ;

[0210] B. The default match score is 2, the default mismatch penalty is 3, the default gap opening penalty is 10, and the default gap extension penalty is 2. Construct a scoring matrix according to the above settings;

[0211] C. Survival cluster identification sequence CS i The reverse complementary (base complementary pairing principle) sequence RCS i ;

[0212] D. Sequence RCS i and reference gene sequence RS i The Smith-Waterman algorithm was used for alignment. If the sequence similarity was greater than SIM (optional parameter, 90% by default), the sequence cluster was considered to be an artificial chimeric sequence cluster of FFPE samples;

[0213] (7) According to the above method, the information set of all artificial chimeric sequence clusters AC is obtained. i ,i∈S *}(S≤M), each artificial chimeric sequence cluster information ACi contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and all artificial chimeric reads {rs j ,j∈K i}.

[0214] 7. The above artificial chimeric sequence clusters are grouped according to the chromosome number chr i And the genome alignment position pos i Merge into a region set AR containing multiple artificial chimeric sequence clusters = {AR k ,k∈T *}, where AR k ={AC i ,i∈R k}, The specific merging process is as follows:

[0215] (1) First, all artificial chimeric sequence clusters are sorted from small to large according to chromosome number and genome alignment position;

[0216] (2) Then, iterate from the first artificial chimeric cluster to the last artificial chimeric cluster. The iterative process is as follows:

[0217] A. Assume that a person's chimeric sequence cluster is AC i , whose chromosome number is chr i , the genome alignment position is pos i ;

[0218] B. The next human chimeric sequence cluster is AC i+1 , whose chromosome number is chr i+1 , the genome alignment position is pos i+1 , if chr i =chr i+1 And pos i +PE≥pos i+1 , then AC i+1 Incorporation into AC i The artificial chimeric cluster region AR k ;

[0219] C.AR k Contains the chromosome number chr k (equal to chr i ), genome alignment starting position start k =min(pos i ,…,pos R ), genome alignment end position end k =max(pos i ,…,pos R ) and the corresponding set of all artificial chimeric sequence clusters {AC i ,i∈R k};

[0220] D. If an artificial chimeric sequence cluster AC j If the previous artificial chimeric sequence cluster does not meet condition B, then AC j As a new artificial chimeric region AR k+1 The next iteration starts with the element; if AR k If the number of artificial chimeric sequence clusters exceeds a certain number (the default is 100, which can be modified), the next artificial chimeric sequence cluster AC will also be included. j As a new artificial chimeric region AR k+1 The next iteration starts with the element.

[0221] 8. The program is designed according to the monitoring mode, and a computing unit is used to monitor the queue information. All the artificial chimeric sequence cluster regions are divided according to the number of processes to form multiple parallel artificial chimeric read segment recognition computing units.

[0222] 9. For each artificial chimeric read segment identification unit, complete the following calculation process:

[0223] (1) Traverse and match to the chromosome number chr k , the genome alignment position is at start k and end k All reads between 1 and 2 are filtered out, and the following types of reads are filtered out: reads that are not aligned to the reference genome, secondary alignments, supplementary alignments, marked as duplicates, or have an alignment quality less than MIN_MQ.

[0224] (2) The matrix M stores the fully matched read base positions (first column) and reference genome positions (second column) by the read cigar value, MD label, and read alignment start position. e ;

[0225] (3) According to the read name, if the read is located in the artificial chimeric cluster library corresponding to the artificial chimeric region, the following judgment is made and the corresponding information is put into the monitoring queue:

[0226] A. If a read is in more than two artificial chimeric clusters, and the clustering direction label ctag of the corresponding artificial chimeric cluster contains both "+" and "-", the read is defined as a non-aligned sequence, and the read name, label rtag and read position of the chimeric site are put into the monitoring queue;

[0227] B. If the ctag of the read corresponding to the artificial sequence cluster is "+", then through the matrix M eGet the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the minimum value as the breakpoint position of the read. Determine the ratio of the length from the position to the 3' end of the read to the entire read length. If it is greater than or equal to a certain value (the default is 0.5), it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number chr k Put it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0228] C. If the ctag of the read corresponding to the artificial sequence cluster is "-", then through the matrix M e Get the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the maximum value as the breakpoint position of the read. Determine the ratio of the length from the 5' end of the read to the position to the entire read length. If it is greater than or equal to a certain value (the default is 0.5), it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number chr k Put it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0229] (4) For the matrix M e , filter the submatrix FM containing the breakpoint information of the artificial chimeric cluster set e If the number of submatrix rows is greater than 0, it is further determined whether the read is an artificial chimeric read. The main step is a multi-threaded judgment process based on the division of artificial chimeric cluster sets:

[0230] A. If the ctag of the artificial chimeric cluster is "+", check whether the length L of the sequence from the read position of the read breakpoint to the 3' end of the read is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence CS i A partial base sequence of length L starting from the 5' end;

[0231] B. If the ctag of the artificial chimeric cluster is "-", check whether the length L of the sequence from the 5' end of the read to the read position of the breakpoint is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence CS i A partial sequence of length L starting from the 3' end of

[0232] If the above conditions are met, the read is considered an artificial chimeric read;

[0233] (5) If after screening, the read meets the artificial chimeric cluster judgment criteria, the information is stored according to the ABC modes of step (3).

[0234] 10. Extract information from the monitoring queue, use the chromosome number as the first-level index, and the read name and read rtag as the second-level index to build an artificial chimeric read segment database. For different data, store the data corresponding to the index as follows:

[0235] For the queue information whose read name is consistent with the read rtag, if there is information marked as non-matching, the read is added with the non-matching tag and stored in the database;

[0236] If the ctag is both "+" and "-", the read will also be tagged with a non-matching tag and stored in the database;

[0237] If only "+" or "-" exists, the queue information is sorted from small to large according to the read position of the chimeric site. If the ctag is "+", the first information is stored. If the ctag is "-", the label ctag of the last queue information and the read position of the chimeric site are stored.

[0238] (II) Identification of false positive mutations in VCF files

[0239] This module builds a process pool according to the artificial sequence cluster recognition units divided in module (I), and uses parallel computing to determine whether a mutation is a false positive mutation caused by an artificial chimeric sequence. If it is caused by an artificial chimeric sequence, the VCF information will be modified and the new information will be stored in the result file. If there is no impact, the original information will be retained in the result file.

[0240] 1. Perform parallel calculations according to the calculation units divided in module (I) to obtain the P value of whether each mutation is affected by the chimeric sequence and various data of mutations excluding the influence of the chimeric sequence:

[0241] (1) Assume that the unit is calculated as unit i, and the chromosome is chr i , the starting position is start i , the end position is end i , including M i mutation sites;

[0242] (2) Chromosome number chr i As an index, query the artificial chimeric read database to store the chimeric read information set corresponding to the chromosome;

[0243] (3) Construct an M i ×10 matrix Matrx i Used to store various types of data for each mutation site, including: the number of reads that support the ref type and are read1, the number of reads that support the ref type and are read2, the number of reads that support the ref type and are read1 and are not artificial chimeric reads, the number of reads that support the ref type and are read2 and are not artificial chimeric reads, the number of reads that support the alt type and are read1, the number of reads that support the alt type and are read2, the number of reads that support the alt type and are read1 and are not artificial chimeric reads, the number of reads that support the alt type and are read2 and are not artificial chimeric reads, the sequencing depth of the site, and the depth of the site excluding artificial chimeric sequences;

[0244] (4) Iterate the chromosome number. The chromosome is chr i , starting position start i and end position end i All reads between the two sequences were filtered out: those that were not aligned to the reference genome, those that were secondary aligned, those that were supplementally aligned, those that were marked as duplicates, or those that had an alignment quality less than MIN_MQ;

[0245] (5) For the retained reads, use the matrix Matrx i Complete multi-threaded judgment, first for each row mutation m (m∈M i ), the mutation information is chr-pos-ref-alt, assuming that the first base of ref and alt is the same, the starting position of the mutation altstart is set to pos-1, if not, altstart is set to pos-2, and the end position of the mutation altend is set to pos plus the length of ref;

[0246] (6) Whether the partial sequence of the read at the mutation position is a chimeric sequence is determined by first using the cigar value and the start position of the read chromosome alignment to obtain the read positions corresponding to altstart and altend: readstart and readend; if the read name plus rtag exists in the artificial chimeric read library and meets the following conditions, it is marked as a chimeric sequence:

[0247] A. Reads are marked as non-aligned;

[0248] B. The ctag of the read is "+", and the read position rpos of the chimeric breakpoint is less than or equal to readstart;

[0249] C. The ctag of Read is "-", and the read position rpos of the chimeric breakpoint is greater than or equal to readend;

[0250] (7)Matrx i The mth row and the 9th column are incremented by 1. If it is judged not to be a chimeric sequence, the mth row and the 10th column are incremented by 1;

[0251] (8) Determine whether the sequence between readstart and readend of read is consistent with ref or alt. If it is consistent with ref, it is marked as isref; if it is consistent with alt, it is marked as isalt.

[0252] (9) If read is read1, set the column position col to 1. If read is read2, set col to 2. If read is marked as isalt, add 4 to the col position. The matrix Matrx i Add 1 to the mth row and colth column;

[0253] (10) If it is not a chimeric sequence, add 2 to col and the matrix Matrx i Add 1 to the mth row and colth column;

[0254] (11) Through the above method, all reads in the calculation unit i area are traversed to obtain the matrix Matrx i ;

[0255] (12) For each mutation, assume that the mutation is located in Matrx i If there are m rows in the matrix, a new two-dimensional matrix Data is constructed, and its calculation rules are as follows:

[0256]

[0257] (13) Data was used as input data and one-sided Fisher's exact test was used to calculate the oddratio (odds ratio) and P value;

[0258] (14) If oddratio is less than 0.5 (default value, which can be modified), the mutation is considered to be affected by the artificial chimeric sequence. Under this condition, if the P value is less than 0.05 (default value, which can be modified), the mutation is considered to be significantly affected by the artificial chimera.

[0259] (15) Construct storage information EAVar for the affected mutations i The storage information is indexed by chr-pos-ref-alt, and the value contains Matrx i [m][3], Matrx i [m][4]、Matrxi [m][7], Matrx i [m][8]、Matrx i [m]

[10] , oddratio, and P value.

[0260] 2. Modify the VCF file, add the ISAE tag to the INFO part of the header file to indicate that the mutation is affected by the chimeric sequence, add the isAEM tag to the FILTER part to indicate that the mutation is a false point caused by the chimeric sequence and is filtered, and then modify the mutation information of each computing unit. Assuming there are K computing units and each computing unit runs in parallel, the modification steps are as follows:

[0261] (1) Prefer to retain more realistic candidate detection mutations, and retain mutations that meet the following conditions: FILTER is PASS, mutation AD (allele depth, mutation support number) ≥ 4 (default value, adjustable), DP (site depth, site sequencing depth) ≥ 40 (default value, adjustable), AF (allele frequency, mutation frequency) ≥ 3% (default value, adjustable);

[0262] (2) Store the total number of mutations Num retained by each computing unit i ;

[0263] (3) Then determine whether the retained mutation is located in EAVar i In the database, if located, the ISAE label is added to the INFO information section to indicate that the mutation is affected by an artificial chimeric sequence;

[0264] (4) Store the total number of mutations INum affected by the artificial chimeric sequence in each computing unit i ;

[0265] (5) According to EAVar i Database information, modify the mutated AD, DP and AF information, and the calculation rules are as follows;

[0266] AD=Matrx i [m][7]+Matrx i [m][8]

[0267] DP=Matrx i [m]

[10]

[0268]

[0269] (6) Then determine whether the mutation meets the conditions in 1). If not, and the P value is less than 0.05, and the FILTER information is PASS, then change the FILTER to the isAEM label;

[0270] (7) Store the total number of mutations FNum that are filtered out due to artificial chimeric sequences in each computing unit i .

[0271] 3. Merge the VCF files of each calculation unit mentioned above to generate a new VCF file, summarizing the mutation information affected by all artificial chimeric sequences.

[0272] 4. Calculate the mutation ratio IPercent affected by the artificial chimeric sequence of the sample and the false point ratio FPercent filtered out. The calculation rules are as follows:

[0273]

[0274]

[0275] IPercent and FPercent are used as quality control indicators to evaluate whether mutations in the next-generation sequencing of FFPE samples are affected by artificial chimeric sequences. When both quality control indicators are greater than 50%, the quality of the FFPE sample is poor.

[0276] Example 2 Comprehensive correction mode for small target area capture

[0277] If sequencing data is captured only for target regions below the Mb level, and the impact of chimeric sequences on sequencing data needs to be completely eliminated, the comprehensive correction mode can be used. The comprehensive correction mode can identify chimeric sequence breakpoints in the comparison data, mark the chimeric sequences as soft clips, and mark the base sequences as N bases to eliminate the impact on downstream analysis.

[0278] (I) Construction of artificial chimeric read database

[0279] This module involves the construction of an artificial chimeric read database. This module is mainly divided into two parallel computing parts: artificial chimeric sequence cluster region identification and artificial chimeric read database construction. This module contains the following contents:

[0280] 1. FFPE samples numbered T3 and T4 undergo nucleic acid extraction, enzyme shearing, library preparation, hybridization capture and sequencing to generate sequencing data. Hybridization capture can use a combination of gene probes for the researcher's interest;

[0281] 2. Use bwa to align the sequencing data with the hg19 reference genome to obtain an alignment file. Use picard to mark the repeated fragments in the alignment file to obtain an alignment file with marked repeated fragments.

[0282] 3. Divide the potential parallel computing units according to different data or computing modes, and capture the sequencing probe combination sequencing for the target area. All target areas can be used as potential parallel computing units.

[0283] 4. Combine the above potential parallel computing units in the following way to reasonably reduce the overall running time:

[0284] (1) Assume that a potential parallel computing unit i, chromosome is chr i , the starting position is start i , the end position is end i , another potential parallel computing unit j, chromosome is chr j , the starting position is start j , the end position is end j , if the above two units meet any of the following conditions, the two calculation units are merged:

[0285] A.chr i =chr j And start i -end j ≤PE;

[0286] B.chr i =chr j And start j -end i ≤PE;

[0287] The above PE is the read length of the paired-end sequencing. For example, if the paired-end sequencing is 150bp, the PE is 150.

[0288] (2) If the number of potential computing units (Num) corresponding to the merged computing unit exceeds a certain value, such as 50, the computing units are split according to the coordinate position to ensure that Num has a certain upper limit and the running time of each merged computing unit is relatively balanced;

[0289] 5. For each computing unit, parallel processing is performed to find all artificial chimeric sequence clusters in the region. The main steps are as follows:

[0290] (1) Traverse all reads in the computational unit region and filter out reads that are not aligned to the reference genome, secondary aligned, supplementary aligned, marked as duplicates, or have an alignment quality less than MIN_MQ (the default value is 20, which can be adjusted);

[0291] (2) The following three matrices are obtained by using the read alignment starting position, MD (a set of strings that record mismatched bases and reference base missing) labels, and cigar (a description of the pattern of sequence alignment to the reference sequence and the number of bases corresponding to each type of pattern) value:

[0292] A. Matrix M storing the base positions of the reads that match perfectly (first column) and the reference genome positions (second column) e ;

[0293] B. Matrix M storing mismatched read base positions (first column) and reference genome positions (second column) s ;

[0294] C. Matrix M that stores the cigar pattern number (first column) and the number of bases corresponding to the cigar pattern (second column) c .

[0295] (3) Artificial chimeric read candidate sequences show the following six patterns in the Bam (Binary Alignment Map) file:

[0296] A. Mode 1: There are more than CLIPLEN (optional parameter, default is 6) continuous softclip bases at the 5' end of the read segment (DNA polymerase working direction, from 5' to 3' end);

[0297] B. Mode 2: There are more than CLIPLEN consecutive softclip bases at the 3' end of the read segment;

[0298] C. Mode 3: There are more than CLIPLEN continuous softclip bases at the 5' end of the read segment, and there is a SNV or INDEL mutation within the downstream ALTDIS (optional parameter, default is 30) bases;

[0299] D. Mode 4: There are more than CLIPLEN continuous softclip bases at the 3' end of the read segment, and there is a SNV or INDEL mutation within ALTDIS bases upstream;

[0300] E. Mode 5: There are more than MINSNV (optional parameter, default is 3) SNV or INDEL mutations within the 5' end of the read segment ALTDIS;

[0301] F. Mode 6: There are MINSNV or more SNV or INDEL mutations within the bases within ALTDIS at the 3' end of the read segment.

[0302] (4) For the candidate chimeric reads of the six modes in (3), the device records the following information through the three matrices in (2):

[0303] A. The reference genome chromosome number chr to which the reads are aligned;

[0304] B. Obtain the read position rpos and the genome alignment position pos of the read chimeric breakpoint. The corresponding processing steps are as follows:

[0305] a. If it is mode 1, the read position of the base 1 bp downstream of the softclip sequence and the corresponding genome alignment position are recorded;

[0306] b. If it is mode 2, the read position of the base 1 bp upstream of the softclip sequence and the corresponding genome alignment position are recorded;

[0307] c. If it is mode 3, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases downstream of the 5' softclip sequence of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0308] d. If it is mode 4, obtain the farthest mutation base (SNV or INDEL) within ALTDIS bases upstream of the 3' end softclip sequence of the read segment B i , search upstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of

[0309] e. If it is mode 5, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases downstream of the 5' end of the read segment B i , search downstream from this position to obtain the first read base B that matches the reference base j If the first base of the read segment is at B j The number of bases between the read positions is greater than CLIPLEN, then record B j The read positions and genome comparison positions of

[0310] f. If it is mode 6, obtain the farthest mutant base (SNV or INDEL) within ALTDIS bases upstream of the 3' end of the read segment B i , search upstream from this position to obtain the first read base B that matches the reference base jIf you read segment B j The number of bases between the read position and the end of the read is greater than CLIPLEN, record B j The read positions and genome comparison positions of

[0311] C. The name of the read segment readname and the tag rtag indicating whether the read segment is Read1 (sequencing sequence of the first round of template strand in double-end sequencing) or Read2 (sequencing sequence of the second round of template strand in double-end sequencing);

[0312] D. Clustering direction tag ctag, for modes 1, 3 and 5, record the tag “-”, indicating that the subsequent clustering direction is from downstream to upstream; for modes 2, 4 and 6, record the tag “+”, indicating that the subsequent clustering direction is from upstream to downstream;

[0313] E. Chimeric sequence S, if the clustering label is "-", the sequence from the 5' end of the read to rpos is recorded; if it is "+", the sequence from rpos to the 3' end of the read is recorded.

[0314] (5) Identifying candidate chimeric read clusters through clustering algorithms:

[0315] A. Prefer to sort by chromosome chr, genome alignment position pos and cluster label;

[0316] B. Cluster the candidate chimeric reads that have the same chromosome chr, genome alignment position pos, and clustering tag ctag. Suppose the set R = {r i ,i∈N *}, r i Recording the candidate chimeric sequence information;

[0317] C. The set R is sorted from large to small according to the length of the chimeric sequence S. Assume that the sorted set is RS = {rs i ,i∈N *};

[0318] D. Starting from rs1 to rs N Until then, iteratively calculate all candidate chimeric sequence cluster sets C, C = {C i ,i∈M *}, where C i = {rs j ,j∈K i}, The specific calculation process is as follows:

[0319] a. Assume that the candidate chimeric read of the starting iteration is rs j , the chimeric sequence is S j , the length of the chimeric sequence is L j , Sj It is a set C i Identification sequence CS of candidate chimeric sequence cluster i ;

[0320] b. From rs j+1 Start, rs j+1 The chimeric sequence is S j+1 , the length of the chimeric sequence is L j+1 , if any of the following conditions are met: If the clustering direction label ctag is "-", S j+1 Equal to S j L j -L j+1 +1 to L j If ctag is "+", S j+1 Equal to S j 1st to Lth j+1 The base sequence between rs j+1 Included in set C j ;

[0321] c. If the kth rs k If the above conditions are not met, rs k The chimeric sequence S k As a new round of cluster set C i+1 The recognition sequence CS i+1 , start the next round of iteration;

[0322] E. Finally, we get a set C, where each element C i Contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and which contains all candidate chimeric reads {rs j ,j∈Ki};

[0323] (6) Using local comparison to determine whether the candidate chimeric sequence cluster is a chimeric sequence, the main process is as follows:

[0324] A. Obtain the candidate chimeric sequence cluster genome alignment position pos i Upstream and downstream FLANK (optional parameter, default is 500, the sum of the front and back extended sequences is 1000bp, it is generally believed that more than 1000bp may be caused by structural variation) reference gene sequence RS between bases i ;

[0325] B. The default match score is 2, the default mismatch penalty is 3, the default gap opening penalty is 10, and the default gap extension penalty is 2. Construct a scoring matrix according to the above settings;

[0326] C. Survival cluster identification sequence CS i The reverse complementary (base complementary pairing principle) sequence RCS i ;

[0327] D. Sequence RCS i and reference gene sequence RS i The Smith-Waterman algorithm was used for alignment. If the sequence similarity was greater than SIM (optional parameter, 90% by default), the sequence cluster was considered to be an artificial chimeric sequence cluster of FFPE samples;

[0328] (7) According to the above method, the information set of all artificial chimeric sequence clusters AC is obtained. i ,i∈S *}(S≤M), each artificial chimeric sequence cluster information AC i Contains the chromosome number chr i , genome alignment position pos i , clustering direction ctag i , cluster identification sequence CS i and all artificial chimeric reads {rs j ,j∈Ki}.

[0329] 6. The above artificial chimeric sequence clusters are grouped according to the chromosome number chr i And the genome alignment position pos i Merge into a region set AR containing multiple artificial chimeric sequence clusters = {AR k ,k∈T *}, where AR k ={AC i ,i∈R k}, The specific merging process is as follows:

[0330] (1) First, all artificial chimeric sequence clusters are sorted from small to large according to chromosome number and genome alignment position;

[0331] (2) Then, iterate from the first artificial chimeric cluster to the last artificial chimeric cluster. The iterative process is as follows:

[0332] A. Assume that a person's chimeric sequence cluster is AC i , whose chromosome number is chr i , the genome alignment position is pos i ;

[0333] B. The next human chimeric sequence cluster is AC i+1 , whose chromosome number is chr i+1 , the genome alignment position is pos i+1 , if chr i =chr i+1 And pos i +PE≥pos i+1 , then AC i+1 Incorporation into AC i The artificial chimeric cluster region AR k ;

[0334] C.AR k Contains the chromosome number chr k (equal to chr i ), genome alignment starting position start k =min(pos i ,…,pos R ), genome alignment end position endk = max(pos i ,…,pos R ) and the corresponding set of all artificial chimeric sequence clusters {AC i ,i∈R k};

[0335] D. If an artificial chimeric sequence cluster AC j If the previous artificial chimeric sequence cluster does not meet condition B, then AC j As a new artificial chimeric region AR k+1 The next iteration starts with the element; if AR k If the number of artificial chimeric sequence clusters exceeds a certain number (the default is 100, which can be modified), the next artificial chimeric sequence cluster AC will also be included. j As a new artificial chimeric region AR k+1 The next iteration starts with the element.

[0336] 7. The program is designed according to the monitoring mode, and a computing unit is used to monitor the queue information. All the artificial chimeric sequence cluster regions are divided according to the number of processes to form multiple parallel artificial chimeric read segment recognition computing units.

[0337] 8. For each artificial chimeric read segment identification unit, complete the following calculation process:

[0338] (1) Traverse and match to the chromosome number chr k , the genome alignment position is at start k and end kAll reads between 1 and 2 are filtered out, and the following types of reads are filtered out: reads that are not aligned to the reference genome, secondary alignments, supplementary alignments, marked as duplicates, or have an alignment quality less than MIN_MQ.

[0339] (2) The matrix M stores the fully matched read base positions (first column) and reference genome positions (second column) by the read cigar value, MD label, and read alignment start position. e ;

[0340] (3) According to the read name, if the read is located in the artificial chimeric cluster library corresponding to the artificial chimeric region, the following judgment is made and the corresponding information is put into the monitoring queue:

[0341] A. If a read is in more than two artificial chimeric clusters, and the clustering direction label ctag of the corresponding artificial chimeric cluster contains both "+" and "-", the read is defined as a non-aligned sequence, and the read name, label rtag and read position of the chimeric site are put into the monitoring queue;

[0342] B. If the ctag of the read corresponding to the artificial sequence cluster is "+", then through the matrix M e Get the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the minimum value as the breakpoint position of the read. Determine the ratio of the length from the position to the 3' end of the read to the entire read length. If it is greater than or equal to a certain value (the default is 0.5), it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number chr k Put it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0343] C. If the ctag of the read corresponding to the artificial sequence cluster is "-", then through the matrix M e Get the read position corresponding to the breakpoint genomic position of the artificial sequence cluster. If there are multiple breakpoint positions, select the maximum value as the breakpoint position of the read. Determine the ratio of the length from the 5' end of the read to the position to the entire read length. If it is greater than or equal to a certain value (the default is 0.5), it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number chr kPut it into the monitoring queue, otherwise the read name, label rtag, label ctag, read position of the chimeric site and chromosome number chr k Put it into the monitoring queue;

[0344] (4) For the matrix M e , filter the submatrix FM containing the breakpoint information of the artificial chimeric cluster set e If the number of submatrix rows is greater than 0, it is further determined whether the read is an artificial chimeric read. The main step is a multi-threaded judgment process based on the division of artificial chimeric cluster sets:

[0345] A. If the ctag of the artificial chimeric cluster is "+", check whether the length L of the sequence from the read position of the read breakpoint to the 3' end of the read is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence CS i A partial base sequence of length L starting from the 5' end;

[0346] B. If the ctag of the artificial chimeric cluster is "-", check whether the length L of the sequence from the 5' end of the read to the read position of the breakpoint is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence CS i A partial sequence of length L starting from the 3' end of

[0347] If the above conditions are met, the read is considered an artificial chimeric read;

[0348] (5) If after screening, the read meets the artificial chimeric cluster judgment criteria, the information is stored according to the ABC modes of step (3).

[0349] 9. Extract information from the monitoring queue, use the chromosome number as the first-level index, and the read name and read rtag as the second-level index to build an artificial chimeric read database. For different data, store the data corresponding to the index as follows:

[0350] For the queue information whose read name is consistent with the read rtag, if there is information marked as non-matching, the read is added with the non-matching tag and stored in the database;

[0351] If the ctag is both "+" and "-", the read will also be tagged with a non-matching tag and stored in the database;

[0352] If only "+" or "-" exists, the queue information is sorted from small to large according to the read position of the chimeric site. If the ctag is "+", the first information is stored. If the ctag is "-", the label ctag of the last queue information and the read position of the chimeric site are stored.

[0353] (II) Correction of artificial chimeric reads in the alignment file

[0354] The artificial chimeric read correction module of the alignment file mainly uses the artificial chimeric read database to modify the artificial chimeric sequence information in the alignment file, change its base information to N base, and modify the cigar value and label information to eliminate the influence of artificial chimeric sequences on downstream analysis.

[0355] 1. First, extract the unaligned sequences into a file named unmap.bam using samtools.

[0356] 2. Divide the artificial chimeric read database in (I) into computing units according to a type of index, i.e., chromosome number. Assuming there are N computing units, process the alignment files in parallel:

[0357] (1) Assume that the chromosome number corresponding to a computational unit is chr i , the output comparison file is Part i .bam;

[0358] (2) In Part i The parameter information of this processing is added to the header file information of .bam for recording;

[0359] (3) Traverse and compare to chr i All reads in Part 1, if the read has secondary alignment, supplementary alignment, and duplicate tags, these reads are written to Part 2. i .bam file and skip the following processing steps;

[0360] (4) If the alignment quality of a read is greater than or equal to MIN_MQ (the default is 20, which can be set), the number of such reads is counted, which is ReadNum i ;

[0361] (5) If the paired read of the read is in the artificial chimeric read database and is marked as non-aligned, the read is labeled with the paired read non-aligned label and then written into Part i .bam file;

[0362] (6) If the read is in the artificial chimeric read database and is marked as non-aligned, the read is labeled as non-aligned, the alignment quality is set to 0, the cigar string is set to ".", and then Part i .bam file;

[0363] (7) If the read is in the artificial chimeric read database and is not marked as non-aligned, if ctag is a "+" tag, change the read sequence from the chimeric breakpoint to the end base to N bases, and change the string corresponding to the bases from the chimeric breakpoint to the end in the cigar value to the number of these bases plus the "S" (softclip) tag, and then write it to Part i .bam file;

[0364] (8) If the read is in the artificial chimeric read database and is not marked as non-aligned, if ctag is the "-" tag, the bases from the beginning of the read to the chimeric breakpoint are changed to N bases, and the reference genome alignment start position of the read is set to the reference genome position corresponding to the next base of the chimeric breakpoint. At the same time, the cigar value from the beginning to the chimeric breakpoint base corresponding to the string is changed to the number of these bases plus the "S" (softclip) tag, and then written to Part i .bam file;

[0365] (9) The number of reads in (6), (7), and (8) is AReadNum i ;

[0366] 3. Calculate the percentage of chimeric reads in the sample, APercent, as an indicator for sample quality control evaluation. The calculation method is as follows:

[0367]

[0368] 4. Merge unmap.bam file and all Parts i .bam file to the same output comparison file out.bam, and set the header file of the output file to Part i The header file of the .bam file, use samatools sort to sort the output comparison file out.bam, and use samtools index to build the index file of the output file;

[0369] 5. The sorted output alignment file and index file are used in subsequent downstream analysis, which removes the influence of chimeric sequences.

[0370] Experimental example

[0371] 1. Samples T1 and T2 are pleural effusion and FFPE samples, respectively. They are interrupted by enzyme cutting, and then somatic cell variation detection is performed by conventional methods. 16 mutation sites are randomly selected, and mPCR (Multiplex Polymerase Chain Reaction) primers are designed. Ultrasonic interruption plus mPCR method is used to construct files, and then target region mutation detection is performed. At the same time, the data of enzyme cutting interruption of T1 and T2 are processed according to the method of Example 1 above to exclude false positive mutations caused by artificial chimeric sequences. The mutation sites detected by the mPCR method are taken as positive sites (P), and those not detected are taken as negative sites (N). Comparing the mutations after the treatment of the present invention, the mutations retained by the present invention in P are TP, and the mutations retained by the present invention in N are FP. TPR (true positive rate) and FDR (false positive rate) are statistically analyzed:

[0372]

[0373]

[0374] The statistical results are as follows:

[0375] Sample No. P TP TPR% N FP FDR% T1 16 16 100 0 NA NA T2 0 NA NA 15 0 0

[0376] There was no false positive mutation in the T1 sample, and most of the T2 samples were false positive sites except for one site where the primer design failed. After being processed by the method of the present invention, the TPR of T1 was 100%, and the FDR of T2 samples was 0%.

[0377] 2. Samples T3 and T4 are FFPE samples. They were processed according to the above method. 32 sites were randomly selected to design primers. The method of Example 2 of the present invention was used to generate a comparison file and detect mutations. The ultrasound interruption mPCR detection was used as the standard. The statistical results are as follows:

[0378] Sample No. P TPR% N FDR% T3 31 100 0 NA T4 4 100 28 0

[0379] T3 has a site with low multiple coverage. It can be seen that the proportion of false positive mutations is inconsistent in different FFPE samples. T3 is all real mutations, and the TPR of the present invention is 100%. T4 has 4 real mutations and 28 false positive mutations. The TPR and FDR of the method of the present invention are 100% and 0%, respectively, which proves that the present invention can effectively eliminate the influence of chimeric reads on mutations.

[0380] 3. For T1 and T2, their sample quality evaluation indicators IPercent and FPercent are as follows:

[0381] Sample No. IPercent FPercent T1 3.68 3.16 T2 68.48 58.73

[0382] From the above indicators, it can be seen that the quality of pleural and ascites samples T1 is better, and the quality of FFPE samples T2 is poorer. Most mutation sites are greatly affected by chimeric reads, and the sample quality control fails.

[0383] 4. For T3 and T4, their sample quality assessment indicators APercent are as follows:

[0384] Sample No. APercent T3 4.99 T4 9.42

[0385] From the above indicators, it can be seen that the proportion of chimeric reads in the target region of T3 is small, and the proportion of chimeric reads in the target region of T4 is large, and the sample quality of the two FFPE samples is different.

[0386] Obviously, the above embodiments are merely examples for the purpose of clear explanation, and are not intended to limit the implementation methods. For those skilled in the art, other different forms of changes or modifications can be made based on the above description. It is not necessary and impossible to list all the implementation methods here. The obvious changes or modifications derived therefrom are still within the scope of protection of the invention.

Claims

1. A computing device for identifying artificial chimeric reads in next-generation sequencing data of FFPE samples, characterized in that: The computing device includes three modules, namely, a module for constructing an artificial chimeric read database, a module for identifying false positive mutations in a VCF file, and a module for correcting artificial chimeric reads in an alignment file; The computational device can be used to identify artificial chimeric reads in FFPE samples in two modes: fast computational mode and comprehensive correction mode. The fast calculation mode includes a module for building an artificial chimeric read database and a module for identifying false positive mutations in VCF files. The comprehensive correction mode includes a construction module for the artificial chimeric read database and a correction module for artificial chimeric reads in the alignment file; Different modes are selected for different gene detection ranges. For small probe combinations, the comprehensive correction mode is selected, and for whole exome sequencing or large probe combinations, the fast calculation mode is selected. The artificial chimeric read database construction module is divided into two parallel computing parts: artificial chimeric sequence cluster region identification and artificial chimeric read database construction, including the following contents: S11, sequencing data is generated from FFPE samples after nucleic acid extraction, enzyme shearing, library preparation, hybridization capture and sequencing; S12, the sequencing data is aligned with the hg19 reference genome by bwa to obtain an alignment file, and the repeated fragments in the alignment file are marked by picard to obtain an alignment file with marked repeated fragments; S13, divide the potential parallel computing units according to different computing modes: S131, if it is a fast computing mode, divide the potential parallel computing units according to the mutations retained in the VCF, and the specific operations are as follows: (1) Retain relatively realistic candidate detection mutations, and retain mutations that meet the following conditions: the filter information FILTER is PASS, the mutation fragment support number AD is greater than or equal to a certain value, the site depth DP is greater than or equal to a certain value, and the mutation frequency AF is greater than or equal to a certain value; (2) The region where the above-mentioned retained mutation is located is used as a potential parallel computing unit. Assuming that the chromosome number of the mutation is chr, the mutation position is pos, the reference base is ref, and the mutation base is alt, then the chromosome number of the potential parallel computing unit is chr, the starting position number is pos-1, and the ending position number is pos-1 plus the length of ref. The mutation information is stored in the structure of chr-pos-ref-alt; S132, if it is a comprehensive correction mode, capture the sequencing probe combination sequencing for the target region, and use all the target regions as potential parallel computing units; S14, merging the above potential parallel computing units according to the following method to reasonably reduce the overall running time: (1) Assume that for a potential parallel computing unit i, the chromosome is , the starting position is , the end position is , another potential parallel computing unit j, chromosome is , the starting position is , the end position is , if the above two units meet any of the following conditions, the two calculation units are merged: a. and ; b. and ; The above PE is the read length of paired-end sequencing; (2) If the number of potential parallel computing units Num corresponding to the merged computing unit exceeds a certain value, the computing units are split according to the coordinate position to ensure that Num has a certain upper limit and the running time of each merged computing unit is relatively balanced. If it is a fast computing mode, for the parallel computing unit, the mutation information corresponding to the unit is stored at the same time; S15, parallel processing for each computational unit to find all artificial chimeric sequence clusters in the region; S16, the artificial chimeric sequence clusters are numbered according to chromosome and genome alignment positions Merge into a collection of regions containing multiple artificial chimeric sequence clusters , Among them, k represents the region after the kth artificial chimeric sequence cluster is merged, Represents the total number of all artificial chimeric sequence cluster regions in the FFPE sample, each artificial chimeric sequence cluster region , represents an artificial chimeric sequence cluster in the region, storing sequence-related information, p represents the pth artificial chimeric sequence cluster, R k The number of artificial chimeric sequence clusters contained in the corresponding k-th region, Total number of artificial chimeric sequence clusters ; S17, program design is performed according to the monitoring mode, a computing unit is specifically used for monitoring queue information, and all the artificial chimeric sequence cluster regions are divided according to the number of processes to form multiple parallel artificial chimeric read segment identification computing units; S18, for each artificial chimeric read segment identification calculation unit, complete the following calculation process: S181, traverse and match to the chromosome number is , the genome alignment position is located at and All reads between n and n are filtered out: reads that are not aligned to the reference genome, secondary aligned, supplementary aligned, marked as duplicates, or have an alignment quality less than n1, where n1 is MIN_MQ; S182, a matrix storing the fully matched read base positions and reference genome positions using the read cigar value, MD tag, and read alignment start position ; S183, judging according to the read segment name, if the read segment is located in the artificial chimeric cluster library corresponding to the artificial chimeric sequence cluster region, the following judgment is performed and the corresponding information is put into the monitoring queue: (1) If a read is in more than two artificial chimeric clusters, and the clustering direction label ctag of the corresponding artificial chimeric cluster contains both "+" and "-", the read is defined as a non-aligned sequence, and the read name, label rtag, and read position of the chimeric site are placed in the monitoring queue. The label rtag indicates whether the read segment is read1 or read2. (2) If the clustering direction label ctag of the read corresponding to the artificial chimeric sequence cluster is "+", then the matrix Get the read position corresponding to the breakpoint genomic position of the artificial chimeric sequence cluster. If there are multiple breakpoint positions, select the minimum value as the breakpoint position of the read, and determine the ratio of the length from the position to the 3' end of the read to the entire read length. If it is greater than or equal to a certain value, it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number are Put it into the monitoring queue, otherwise the read name, label rtag, cluster direction label ctag, read position of the chimeric site and chromosome number Put it into the monitoring queue; (3) If the clustering direction label ctag of the read corresponding to the artificial chimeric sequence cluster is "-", then the matrix Obtain the read position corresponding to the breakpoint genomic position of the artificial chimeric sequence cluster. If there are multiple breakpoint positions, select the maximum value as the breakpoint position of the read, and determine the ratio of the length from the 5' end of the read to the position to the entire read length. If it is greater than or equal to a certain value, it is considered that most of the read sequence is an artificial chimeric sequence, and the read is defined as a non-aligned sequence. At the same time, the read name, label rtag and chromosome number are Put it into the monitoring queue, otherwise the read name, label rtag, cluster direction label ctag, read position of the chimeric site and chromosome number Put it into the monitoring queue; S184, for the matrix , filter the submatrix containing breakpoint information of artificial chimeric cluster sets , if the number of submatrix rows is greater than 0, further determine whether the read is an artificial chimeric read: (1) If the clustering direction tag ctag of the artificial chimeric sequence cluster is "+", check whether the length L of the sequence from the read position of the read breakpoint to the 3' end of the read is greater than or equal to n2, where n2 is CLIPLEN, and whether the sequence is equal to the cluster identification sequence A partial base sequence of length L starting from the 5' end; (2) If the clustering direction tag ctag of the artificial chimeric cluster is "-", check whether the length L of the sequence from the 5' end of the read to the read position of the breakpoint is greater than or equal to CLIPLEN and the sequence is equal to the cluster identification sequence A partial sequence of length L starting from the 3' end of If the above conditions are met, the read is considered an artificial chimeric read; S185, if after screening, the read meets the artificial chimeric cluster judgment conditions, then the information is stored according to the three modes (1), (2) and (3) of step S183; S19, extract information from the monitoring queue, use the chromosome number as the first-level index, and the read name and read rtag as the second-level index to build an artificial chimeric read database; for different data, store the data corresponding to the index according to the following method: For the queue information whose read name is consistent with the read rtag, if there is information marked as non-matching, the read is added with the non-matching tag and stored in the database; If the clustering direction tag ctag is both "+" and "-", the read will also be added with a non-matching tag and stored in the database; If only "+" or "-" exists, the queue information is sorted from small to large according to the read position of the chimeric site. If the clustering direction label ctag is "+", the first information is stored. If the clustering direction label ctag is "-", the clustering direction label ctag and the read position of the chimeric site of the last queue information are stored. The false positive mutation recognition module in the VCF file builds a process pool according to the artificial chimeric sequence cluster recognition calculation units divided in the construction module of the artificial chimeric read database, and uses parallel computing to determine whether a mutation is a false positive mutation caused by an artificial chimeric sequence; If it is caused by artificial chimeric sequences, the VCF information will be modified and the new information will be stored in the result file. If there is no impact, the original information will be retained in the result file; The specific steps include: S21, performing parallel calculations according to the calculation units divided in the construction module of the artificial chimeric read database to obtain a P value of whether each mutation is affected by the chimeric sequence and various data of mutations excluding the influence of the chimeric sequence: (1) Assume that the computational unit q, the chromosome is , the starting position is , the end position is ,Include mutation sites; (2) Chromosome numbering As an index, query the artificial chimeric read database to store the chimeric read information set corresponding to the chromosome; (3) Build a The matrix Used to store various types of data for each mutation site, including: the number of reads that support the ref type and are read1, the number of reads that support the ref type and are read2, the number of reads that support the ref type and are read1 and are not artificial chimeric reads, the number of reads that support the ref type and are read2 and are not artificial chimeric reads, the number of reads that support the alt type and are read1, the number of reads that support the alt type and are read2, the number of reads that support the alt type and are read1 and are not artificial chimeric reads, the number of reads that support the alt type and are read2 and are not artificial chimeric reads, the sequencing depth of the site, and the depth of the site excluding artificial chimeric sequences; (4) Iteration chromosome number The chromosome is , starting position and end position All reads between the two sequences were filtered out: those that were not aligned to the reference genome, those that were secondary aligned, those that were supplementally aligned, those that were marked as duplicates, or those that had an alignment quality less than MIN_MQ; (5) For the retained reads, use the matrix Complete multi-threaded judgment, first for each row mutation m( ), the mutation information is chr-pos-ref-alt, assuming that the first base of ref and alt is the same, the starting position of the mutation altstart is set to pos-1, if not, altstart is set to pos-2, and the end position of the mutation altend is set to pos plus the length of ref; (6) Whether the partial sequence of the read at the mutation position is a chimeric sequence, the read positions corresponding to altstart and altend are obtained by using the cigar value and the start position of the read chromosome alignment: readstart and readend; if the read name plus rtag exists in the artificial chimeric read database and meets the following conditions, it is marked as a chimeric sequence: a. Reads are marked as non-aligned; b. The clustering direction tag ctag of the read is "+", and the read position rpos of the chimeric breakpoint is less than or equal to readstart; c. The clustering direction label ctag of the read is "-", and the read position rpos of the chimeric breakpoint is greater than or equal to readend; (7) The mth row and the 9th column are incremented by 1. If it is judged not to be a chimeric sequence, the mth row and the 10th column are incremented by 1; (8) Determine whether the sequence between readstart and readend of read is consistent with ref or alt. If it is consistent with ref, mark it as isref; if it is consistent with alt, mark it as isalt. (9) If read is read1, set the column position col to 1. If read is read2, set col to 2. If read is marked as isalt, add 4 to the col position. The matrix Add 1 to the mth row and colth column; (10) If it is not a chimeric sequence, add 2 to col and the matrix Add 1 to the mth row and colth column; (11) Through the above method, all reads in the calculation unit q area are traversed to obtain the matrix ; (12) For each mutation, assume that the mutation is located m rows, then construct a new two-dimensional matrix , and its calculation rules are as follows: ； (13) As input data, one-sided Fisher's exact test was used to calculate the oddratio and P value; (14) If oddratio is less than a certain value, the mutation is considered to be affected by the artificial chimeric sequence. Under this condition, if the P value is less than a certain value, the mutation is considered to be significantly affected by the artificial chimera. (15) Construct storage information for affected mutations The storage information is indexed by chr-pos-ref-alt and the value contains , , , , , oddratio and P value; S22, modify the VCF file, add the ISAE tag to the INFO part of the header file, indicate that the mutation is affected by the chimeric sequence, add the isAEM tag to the FILTER part, indicate that the mutation is a false point caused by the chimeric sequence and is filtered, and then modify the mutation information of each computing unit. Assuming there are K computing units and each computing unit runs in parallel, the modification steps are as follows: (1) Retain relatively realistic candidate detection mutations, and retain mutations that meet the following conditions: FILTER is PASS, the AD of the mutation is greater than or equal to a certain value, the DP is greater than or equal to a certain value, and the AF is greater than or equal to a certain value; (2) Store the total number of mutations retained for each computational unit ; (3) Determine whether the retained mutation is located in In the database, if located, the ISAE label is added to the INFO information section to indicate that the mutation is affected by an artificial chimeric sequence; (4) Store the total number of mutations affected by the artificial chimeric sequence in each computational unit ; (5) Based on Database information, modify the mutated AD, DP and AF information, and the calculation rules are as follows; (6) Then determine whether the mutation meets the conditions in (1). If it does not meet the conditions, and the P value is less than 0.05, the FILTER information is PASS, then change the FILTER to the isAEM label; (7) Store the total number of mutations filtered out due to artificial chimeric sequences in each computational unit ; S23, merging the VCF files of each of the above calculation units to generate a new VCF file, summarizing the mutation information affected by all artificial chimeric sequences; S24, calculate the mutation ratio affected by the artificial chimeric sequence of the sample And the proportion of false points filtered out , the calculation rules are as follows: Will and As quality control indicators for evaluating whether mutations in next-generation sequencing of FFPE samples are affected by artificial chimeric sequences, when both quality control indicators are greater than 50%, the quality of the FFPE sample is poor; The correction module for artificial chimeric reads in the alignment file uses an artificial chimeric read database to modify the artificial chimeric sequence information in the alignment file, change its base information to N bases, and modify the cigar value and label information to eliminate the influence of the artificial chimeric sequence on downstream analysis, including the following steps: S31, extract the unmapped read information into a file named unmap.bam through samtoolsview; S32, dividing the artificial chimeric read database in the construction module of the artificial chimeric read database into computing units according to the first-level index, i.e., the chromosome number, assuming that there are N computing units, and processing the alignment files in parallel: (1) Assume that the chromosome number corresponding to a computational unit is , the output comparison file is ; (2) In Add the parameter information of this processing to the header file information for recording; (3) Traverse and compare to If the read has a secondary alignment, a supplementary alignment, and a duplicate tag, these reads are written file and skip the following processing steps; (4) If the alignment quality of the read is greater than or equal to MIN_MQ, the number of such reads is counted, which is ; (5) If the paired read of the read is in the artificial chimeric read database and is marked as non-aligned, the read is labeled with the paired read non-aligned label and then written document; (6) If the read is in the artificial chimeric read database and is marked as non-aligned, add the non-aligned label to the read, set the alignment quality to 0, set the cigar string to ".", and then write document; (7) If the read is in the artificial chimeric read database and is not marked as non-aligned, if the clustering direction label ctag is a "+" label, then change the read sequence from the chimeric breakpoint to the end base to N bases, and change the string corresponding to the bases from the chimeric breakpoint to the end in the cigar value to the number of bases plus the label, and then write document; (8) If the read is located in the artificial chimeric read database and is not marked as non-aligned, if the cluster direction label ctag is a "-" label, the bases from the beginning of the read to the chimeric breakpoint are changed to N bases, and the reference genome alignment start position of the read is set to the reference genome position corresponding to the next base of the chimeric breakpoint. At the same time, the cigar value from the beginning to the chimeric breakpoint base corresponding to the string is changed to the number of bases plus the label, and then written document; (9) The number of reads in cases (6), (7), and (8) is ; S33. Calculate the percentage of chimeric reads in this sample , as an indicator for sample quality control evaluation, is calculated as follows: ； S34, merge unmap.bam file and all The files are transferred to the same output comparison file out.bam, and the header file of the output file is set to The header file of the file, use samatoolssort to sort the output comparison file out.bam, and use samtoolsindex to build the index file of the output file; S35. The sorted output alignment file and index file are used in subsequent downstream analysis, and the file removes the influence of chimeric sequences.

2. The computing device according to claim 1, characterized in that: The specific steps of S15 are as follows: S151, traverse all reads in the computational unit region, and filter out reads that are not aligned to the reference genome, are secondary aligned, are supplementally aligned, are marked as duplicates, or have an alignment quality less than MIN_MQ; S152, the following three matrices are obtained by read alignment starting position, MD label and cigar value: (1) A matrix storing the base positions of fully matched reads and reference genome positions ; (2) A matrix storing mismatched read base positions and reference genome positions ; (3) A matrix storing the cigar pattern number and the number of bases corresponding to the cigar pattern ; S153, artificial chimeric read candidate sequences show the following 6 patterns in the Bam file: (1) Mode 1: There are more than CLIPLEN consecutive softclip bases at the 5' end of the read segment; (2) Mode 2: There are more than CLIPLEN consecutive softclip bases at the 3' end of the read segment; (3) Mode 3: There are more than CLIPLEN continuous softclip bases at the 5' end of the read segment, and there is a SNV or INDEL mutation within n3 bases downstream, and n3 is ALTDIS; (4) Mode 4: There are more than CLIPLEN continuous softclip bases at the 3' end of the read segment, and there is a SNV or INDEL mutation within ALTDIS bases upstream; (5) Mode 5: There are more than n4 SNVs or INDEL mutations within the 5' end of the read segment within ALTDIS, n4 is MINSNV, and MINSNV is greater than or equal to 2; (6) Mode 6: There are more than MINSNV SNV or INDEL mutations within the ALTDIS bases at the 3' end of the read segment; S154, for the candidate chimeric reads of the six modes in S153, the following information is recorded through the three matrices in S152: (1) The chromosome number of the reference genome to which the read segment is aligned, chr; (2) Obtain the read position rpos and genome alignment position pos of the read chimeric breakpoint. The corresponding processing steps are as follows: a. If it is mode 1, the read position of the base 1 bp downstream of the softclip sequence and the corresponding genome alignment position are recorded; b. If it is mode 2, the read position of the base 1 bp upstream of the softclip sequence and the corresponding genomic alignment position are recorded; c. If it is mode 3, obtain the farthest mutant base B within ALTDIS bases downstream of the 5' end softclip sequence of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of d. If it is mode 4, obtain the farthest mutant base B within ALTDIS bases upstream of the 3' softclip sequence of the read segment i , search upstream from this position to obtain the first read base B that matches the reference base j , record B j The read positions and genome comparison positions of e. If it is mode 5, obtain the farthest mutant base B within ALTDIS bases downstream of the 5' end of the read segment i , search downstream from this position to obtain the first read base B that matches the reference base j If the first base of the read segment is at B j The number of bases between the read positions is greater than CLIPLEN, then record B j The read positions and genome comparison positions of f. If it is mode 6, obtain the farthest mutant base B within ALTDIS bases upstream of the 3' end of the read segment i , search upstream from this position to obtain the first read base B that matches the reference base j If the read base B j The number of bases between the read position and the end of the read is greater than CLIPLEN, record B j The read positions and genome comparison positions of (3) The name of the read segment readname and the label rtag indicating whether the read segment is Read1 or Read2; (4) Clustering direction tag ctag: for modes 1, 3, and 5, the tag "-" is recorded, indicating that the subsequent clustering direction is from downstream to upstream; for modes 2, 4, and 6, the tag "+" is recorded, indicating that the subsequent clustering direction is from upstream to downstream; (5) Chimeric sequence S: if the clustering direction label is "-", the sequence from the 5' end of the read to rpos is recorded; if it is "+", the sequence from rpos to the 3' end of the read is recorded; S155, identifying candidate chimeric read clusters by clustering algorithm: (1) Sort by chromosome chr, genome alignment position pos and clustering direction label; (2) Cluster the candidate chimeric reads that have the same chromosome chr, genome alignment position pos, and clustering direction tag ctag, assuming the set , Record the candidate chimeric read clusters; (3) Collection Sort the chimeric sequences S from large to small in length. Assume that the sorted set is ; (4) From Start to Until then, iteratively calculate all candidate chimeric read cluster sets , ,in , , the specific calculation process is as follows: a. Assume that the candidate chimeric read for the initial iteration is , the chimeric sequence is , the chimeric sequence length is , Is a collection Identification sequences of candidate chimeric read clusters ; b. From start, The chimeric sequence is , the chimeric sequence length is , if any of the following conditions are met: if the clustering direction tag ctag is "-", equal No. to The base sequence between If the clustering direction tag ctag is "+", equal 1st to 2nd The base sequence between Add to collection ; c. If the kth If the above conditions are not met, chimeric sequence As a new round of cluster collection The recognition sequence , start the next round of iteration; (5) Finally, we get the set , each element Contains chromosome number , genome alignment position , cluster direction label , cluster identification sequence and all candidate chimeric reads ; S156, using local comparison to determine whether the candidate chimeric read cluster is a chimeric sequence, the main process is as follows: (1) Obtaining the genome alignment position of candidate chimeric sequence clusters Reference gene sequence between n5 bases upstream and downstream , n5 is FLANK; (2) The default match score is 2, the default mismatch penalty is 3, the default gap opening penalty is 10, and the default gap extension penalty is 2. The scoring matrix is ​​constructed according to the above settings; (3) Generate cluster identification sequence The reverse complement sequence ; (4) Change the sequence and reference gene sequences The Smith–Waterman algorithm was used for alignment, and if the sequence similarity was greater than n6, where n6 was SIM, the sequence cluster was considered to be an artificial chimeric sequence cluster of FFPE samples; S157, obtaining all artificial chimeric sequence cluster information sets according to the above method , each artificial chimeric sequence cluster information Contains chromosome number , genome alignment position , cluster direction label , cluster identification sequence and it contains all artificial chimeric reads .

3. The computing device according to claim 1, characterized in that: In S16, the specific merging process is as follows: S161, sort all artificial chimeric sequence clusters from small to large according to chromosome number and genome alignment position; S162, iterating from the first artificial chimeric sequence cluster to the last artificial chimeric sequence cluster, the iterative process is as follows: (1) Assume that a certain artificial chimeric sequence cluster is , whose chromosome number is , the genome alignment position is ; (2) The next artificial chimeric sequence cluster is , whose chromosome number is , the genome alignment position is ,if and , then Inclusion Artificial chimeric sequence cluster region ; (3) Contains chromosome number , genome alignment starting position , genome alignment termination position And the corresponding set of all artificial chimeric sequence clusters ; (4) If an artificial chimeric sequence cluster If the previous artificial chimeric sequence cluster does not meet condition (2), then As a new artificial chimeric sequence cluster region The next iteration starts with the element of If the number of artificial chimeric sequence clusters contained exceeds a certain number, the next artificial chimeric sequence cluster will also be As a new artificial chimeric sequence cluster region The next iteration starts with the element.

4. A method for identifying artificial chimeric reads in next-generation sequencing data of FFPE samples, characterized by: The sequencing data is generated after FFPE samples undergo nucleic acid extraction, enzyme shearing, library preparation, hybridization capture and on-machine sequencing. The sequencing data is processed using the computing device described in any one of claims 1 to 3, thereby eliminating false positive mutations caused by artificial chimeric sequences and identifying artificial chimeric reads in the second-generation sequencing data of FFPE samples.

Images