Gene typing assembly method and device, electronic equipment and storage medium

Abstract

The embodiment of the invention provides a genotyping assembly method and device, electronic equipment and a storage medium, and belongs to the technical field of biological information. The genotyping assembly method comprises: obtaining target comparison data from sequencing data of a target genome sequence, the target comparison data being used for representing comparison data of any two target fragment sequences in the target genome sequence; a path planning model is constructed according to the target comparison data and the target fragment sequence, and the path planning model is used for representing a path for carrying out typing assembly on the target fragment sequence according to the original sequence data; and taking corresponding original sequence data when the path of the path planning model is minimized as target sequence data, and performing typing assembly on the target fragment sequence according to the target sequence data. According to the embodiment of the invention, genotyping and genomic packaging can be processed at the same time.

Description

Technical Field

[0001] The present disclosure relates to the field of bioinformatics, and in particular to a gene typing assembly method and device, an electronic device and a storage medium. Background Art

[0002] Genotyping is the process of checking DNA sequences and comparing target sequences with reference sequences to determine individual genetic differences. Genome assembly is the process of splicing and assembling the read fragments generated by sequence sequencing to generate the base sequence of the genome. In related technologies, genome sequencing data analysis tools use a distributed processing strategy when processing genotyping and gene assembly, that is, clearly distinguishing the processes of genotyping and assembly. Therefore, there is currently a lack of a method that can handle genotyping and gene assembly at the same time. Summary of the invention

[0003] The embodiments of the present disclosure provide a genotyping assembly method and apparatus, an electronic device, and a storage medium, which are intended to simultaneously process genotyping and gene assembly.

[0004] To achieve the above objectives, the first aspect of the embodiments of the present application proposes a genotyping assembly method, the method comprising:

[0005] Obtaining target comparison data from the sequencing data of the target genome sequence, wherein the target comparison data is used to represent the comparison data of any two target fragment sequences in the target genome sequence;

[0006] Constructing a path planning model according to the target comparison data and the target fragment sequence, wherein the path planning model is used to represent a path for typing and assembling the target fragment sequence according to the original sequence data;

[0007] The original sequence data corresponding to the path minimization of the path planning model is used as the target sequence data, and the target fragment sequence is typed and assembled according to the target sequence data.

[0008] In some embodiments, any two of the target fragment sequences are taken as a target sequence pair;

[0009] The step of obtaining target comparison data from the sequencing data of the target genome sequence comprises:

[0010] Obtaining original comparison data from the sequencing data, wherein the original comparison data is used to represent the comparison data of any two original fragment sequences in the target genome sequence, the original comparison data includes the original comparison starting position and the original overlapping region length, and taking any two original fragment sequences as an original sequence pair;

[0011] Taking any original fragment sequence in the original sequence pair as a comparison sequence, if the original comparison start position is in the front region of the comparison sequence, taking the original sequence pair as a candidate sequence pair;

[0012] Screening the candidate sequence pairs according to the length of the original overlapping region to obtain the target sequence pairs;

[0013] The original alignment data corresponding to the target sequence pair is used as the target alignment data.

[0014] In some embodiments, the original comparison starting position includes a first starting position and a second starting position, and if the original comparison starting position is in the front region of the comparison sequence, taking the original sequence pair as a candidate sequence pair includes:

[0015] Compare the first starting position and the second starting position with a preset position threshold value respectively to obtain a position comparison result;

[0016] If the position comparison result indicates that the first starting position or the second starting position is smaller than the preset position threshold, the original sequence pair is used as a candidate sequence pair.

[0017] In some embodiments, the original alignment data further includes a first sequence length and a second sequence length;

[0018] The step of screening the candidate sequence pairs according to the length of the original overlapping region to obtain the target sequence pairs comprises:

[0019] Calculate the first sequence length by multiples according to a first preset multiple to obtain a first multiple length;

[0020] Calculate the multiple of the second sequence length according to the first preset multiple to obtain a second multiple length;

[0021] Compare the original overlapping area length with the first multiple length and the second multiple length respectively to obtain a first length comparison result;

[0022] If the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, the target sequence pair is obtained according to the candidate sequence pair.

[0023] In some embodiments, if the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, obtaining the target sequence pair according to the candidate sequence pair includes:

[0024] If the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, taking the candidate sequence pair as an initial sequence pair;

[0025] Calculating the first sequence length by multiples according to a second preset multiple to obtain a third multiple length, wherein the second preset multiple is greater than the first preset multiple;

[0026] Calculate the multiple of the second sequence length according to the second preset multiple to obtain a fourth multiple length;

[0027] Compare the original overlapping area length with the third multiple length and the fourth multiple length respectively to obtain a second length comparison result;

[0028] If the second length comparison result indicates that the original overlapping region length is greater than or equal to the third multiple length, and the original overlapping region length is greater than or equal to the fourth multiple length, the initial sequence pair is used as the target sequence pair.

[0029] In some embodiments, constructing a path planning model according to the target alignment data and the target fragment sequence includes:

[0030] Calculate the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights;

[0031] The path planning model is constructed according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence quantity of the target fragment sequence, and the preset haplotype quantity.

[0032] In some embodiments, the target comparison data includes a target overlapping region length and a target number of bases, and the preset weights include a first weight and a second weight;

[0033] The step of calculating the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights includes:

[0034] Perform product calculation based on the first weight and the target base number to obtain first prompt data;

[0035] Calculating the ratio of the target overlapping region length to the target base number to obtain base ratio data;

[0036] Perform product calculation based on the second weight and the base ratio data to obtain second prompt data;

[0037] The typing assembly prompt data is obtained by performing a sum calculation based on the second prompt data, the first prompt data and the target overlapping area length.

[0038] In some embodiments, the target alignment data includes the target overlap region length and the target base number;

[0039] The step of calculating the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights includes:

[0040] Perform product calculation on the target number of bases according to the preset weights to obtain base impact data;

[0041] Performing difference calculation based on the target overlapping region length and the base impact data to obtain third prompt data;

[0042] A fourth prompt data is obtained by performing a sum calculation based on the target overlapping region length and the base influence data;

[0043] The typing and assembly prompt data are obtained by performing a ratio calculation based on the third prompt data and the fourth prompt data.

[0044] In some embodiments, before taking the original sequence data corresponding to the path minimization of the path planning model as the target sequence data, the method further includes performing path minimization processing on the path planning model, specifically including:

[0045] Constructing constraint conditions, wherein the constraint conditions include that the number of typing and assembling of each target fragment sequence is equal to 1, only one target fragment sequence can be typed and assembled at a time, and the typing and assembly starts from a preset start and end sequence;

[0046] The path planning model is subjected to path minimization processing according to the constraint conditions.

[0047] In some embodiments, before constructing the constraint data, the method further includes constructing a preset start and end sequence, specifically including:

[0048] Construct the original start-end sequence;

[0049] The original start-end sequence is compared with the target fragment sequence. If the original start-end sequence and the target fragment sequence do not include an overlapping region, the original start-end sequence is used as the preset start-end sequence.

[0050] In some embodiments, constructing the path planning model according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence number of the target fragment sequence, and the preset number of haplotypes includes:

[0051] Constructing an initial model according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence quantity of the target fragment sequence, and the preset number of haplotypes;

[0052] Determining the chromosome homology status of any two target fragment sequences;

[0053] Constructing the path planning model according to the homologous state of the chromosome and the initial model;

[0054] The constraint condition also includes that the chromosome homology state indicates that any two target fragment sequences originate from the same homologous chromosome.

[0055] In some embodiments, performing path minimization processing on the path planning model according to the constraint condition includes:

[0056] Get the annealing time and annealing times;

[0057] The path planning model is subjected to path minimization processing according to a preset quantum annealing model, the annealing duration, the number of annealing times, and the constraint conditions.

[0058] To achieve the above-mentioned purpose, the second aspect of the embodiment of the present application proposes a genotyping assembly device, the device comprising:

[0059] A comparison data acquisition module, used to obtain target comparison data from the sequencing data of the target genome sequence, wherein the target comparison data is used to represent the comparison data of any two target fragment sequences in the target genome sequence;

[0060] A path planning model construction module, used to construct a path planning model according to the target comparison data and the target fragment sequence, wherein the path planning model is used to represent a path for typing and assembling the target fragment sequence according to the original sequence data;

[0061] The typing and assembly module is used to use the original sequence data corresponding to the path minimization of the path planning model as the target sequence data, and to perform typing and assembly on the target fragment sequence according to the target sequence data.

[0062] To achieve the above-mentioned purpose, the third aspect of an embodiment of the present application proposes an electronic device, comprising a memory and a processor, wherein the memory stores a computer program, and when the processor executes the computer program, the genotyping assembly method described in the first aspect is implemented.

[0063] To achieve the above-mentioned purpose, the fourth aspect of the embodiments of the present application proposes a computer-readable storage medium, which stores a computer program. When the computer program is executed by a processor, it implements the genotyping assembly method described in the first aspect.

[0064] To achieve the above-mentioned purpose, the fifth aspect of the embodiments of the present application proposes a computer program product, which includes a computer program, and the computer program is read and executed by a processor of a computer device, so that the computer device executes the genotyping assembly method described in the first aspect.

[0065] The genotyping assembly method and device, electronic device and storage medium proposed in the embodiment of the present application, build a path planning model through target comparison data and target fragment sequence, and the path planning model can be used to describe the path of typing and assembling the target fragment sequence according to the original sequence data. It can be seen from this that solving the path planning model is actually a combinatorial optimization problem, and the corresponding target sequence data can be used as the optimal solution of the combinatorial optimization problem when the path planning model is minimized. Among them, the path planning model can be processed by path minimization based on the genotyping purpose and the gene assembly purpose. Therefore, when the target fragment sequence is typed and assembled according to the target sequence data, genotyping and gene assembly can be achieved simultaneously. BRIEF DESCRIPTION OF THE DRAWINGS

[0066] The accompanying drawings are used to provide further understanding of the technical solution of the present disclosure and constitute a part of the specification. Together with the embodiments of the present disclosure, they are used to explain the technical solution of the present disclosure and do not constitute a limitation on the technical solution of the present disclosure.

[0067] Figure 1 is a flow chart of the genotyping assembly method provided in the embodiments of the present application;

[0068] Figure 2 yes Figure 1 Flow chart of step S101 in FIG.

[0069] Figure 3 is a schematic diagram showing the original comparison data in the form of a table in the embodiment of the present application;

[0070] Figure 4 is a schematic diagram of the original sequence pairs that need to be filtered out in the embodiments of the present application;

[0071] Figure 5 yes Figure 2 Flow chart of step S202 in FIG.

[0072] Figure 6 yes Figure 2 Flow chart of step S203 in FIG.

[0073] Figure 7 yes Figure 6 Flow chart of step S604 in FIG.

[0074] Figure 8 yes Figure 1 Flow chart of step S102 in FIG.

[0075] Fig. 9 This is a schematic diagram of the comparison of the target fragment sequence position in front and the position in the back in the embodiment of the present application;

[0076] Fig.10 yes Figure 8 Flow chart of step S801 in FIG.

[0077] Fig.11 yes Figure 8 A flowchart of another embodiment of step S801 in FIG.

[0078] Fig.12 is a flowchart of minimizing the path planning model provided by an embodiment of the present application;

[0079] Fig.13 is a flow chart of a method for constructing a preset start and end sequence provided in an embodiment of the present application;

[0080] Fig.14 yes Figure 8 Flowchart of step S802 in FIG.

[0081] Fig.15 yes Fig.12 Flowchart of step S1202 in FIG.

[0082] Fig.16 is a schematic diagram of a simulated triploid short DNA sequence provided in an embodiment of the present application;

[0083] Fig.17 It is a schematic diagram of the sub-type assembly provided in the embodiment of the present application;

[0084] Fig.18 is a schematic diagram of quantum computing operating parameters provided in an embodiment of the present application;

[0085] Fig.19 is a schematic diagram of a genotyping assembly device provided in an embodiment of the present application;

[0086] Fig. 20 It is a schematic diagram of the hardware structure of the electronic device provided in the embodiment of the present application. DETAILED DESCRIPTION

[0087] In order to make the purpose, technical solution and advantages of the present disclosure more clear, the present disclosure is further described in detail below in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain the present disclosure and are not used to limit the present disclosure.

[0088] Before further describing the embodiments of the present disclosure in detail, the nouns and terms involved in the embodiments of the present disclosure are described. The nouns and terms involved in the embodiments of the present disclosure are subject to the following interpretations:

[0089] Genotyping, also known as gene phase analysis, is used to determine the relative order of alleles at multiple loci of an individual, also known as "phase". It indicates which alleles come from the mother and which come from the father. From the perspective of traditional bioinformatics algorithms, typing is to divide sequence fragments (i.e. reads) into two categories and then assemble them separately. The process of genotyping for diploids is a two-cluster problem, and for polyploids it is a multi-cluster problem.

[0090] Genome assembly refers to the use of sequencing methods to generate reads for the genome of the species to be tested, and splicing the fragments according to the overlapping regions between the reads to obtain a higher-quality whole-genome sequence.

[0091] In the related art, genome sequencing data analysis tools adopt a distributed processing strategy when processing genotyping and gene assembly, that is, the processes of genotyping and gene assembly are clearly distinguished. The tools in the related art convert the tasks corresponding to genotyping or gene assembly into complex optimization problems, and the global optimal solution of the optimization problem is the correct genotyping result or the correct gene assembly result. For example, the genotyping task can be converted into minimum error correction (MEC), and the gene assembly can be converted into overlapping layout consensus (OLC).

[0092] Among them, the basic idea of ​​the MEC model is to find a possible combination of haplotypes in a given genotype so that the number of error corrections required when combining these haplotypes is minimized. For diploids, the MEC model is equivalent to the maximum cut problem, which is a typical NP-hard combinatorial optimization problem in graph theory. The points in the graph are reads, and the weighted edges between the points are the number of differences in single nucleotide polymorphism (SNP) sites between reads. By solving the maximum cut problem, the points in the graph are divided into two sets so that the sum of the weights of all edges in the set is minimized, thereby completing the genotyping of diploid genome sequencing data. Commonly used diploid genotyping tools include WhatsHap and HapCUT2. For polyploid genotyping, most tools (such as WhatsHap) cluster reads through clustering algorithms, or tools such as H-PoP and HapTree use the MEC model to solve polyploid genotyping.

[0093] The OLC algorithm starts from reads and reconstructs the original genome sequence by aligning and splicing the reads. The core of the OLC algorithm is to solve the overlapping graph, where the points in the graph are reads, and the weighted edges between the points are the quality of the overlapping areas (overlap) between the reads. A directed path passing through all the points on the graph can be obtained through the optimization problem solving algorithm to determine the relative order of the reads. The process of solving the directed path is equivalent to solving the Hamiltonian path problem in the combinatorial optimization problem.

[0094] In addition, the gene assembly can be converted into the De Bruijn Graph Algorithm (DBG), which also requires the construction of a directed graph. Specifically, the reads are first divided into segments of length k, called k-mers. Each k-mer is represented as a point on the graph, and the edge on the graph represents the connection between two k-mers. On the established directed graph, a path that passes through all edges and only once is found. The found path is the assembly path. The process of finding a path is equivalent to solving the Euler path problem in the combinatorial optimization problem.

[0095] In the related art, heuristic algorithms can be used to solve the above optimization problems. Heuristic algorithms have a short calculation time, but the result quality is low and it is easy to fall into the local optimal solution. Meta-heuristic algorithms can solve the problem of local optimal solutions of heuristic algorithms, but the calculation time of meta-heuristic algorithms is long.

[0096] Quantum computing can solve the problems of heuristic algorithms and metaheuristic algorithms. Quantum computing can process a large amount of information at the same time, improve computing parallelism, and help solve the problem of global optimal solutions. Specifically, quantum optimization methods such as quantum annealing technology and quantum approximate optimization search in the search space through the quantum tunneling effect, which to a certain extent avoids falling into the situation of local optimal solutions. Quantum computing solves the problem of limitations of heuristic algorithms and metaheuristic algorithms in related technologies. Quantum computing can not only find the global optimal solution, but also avoid the problem of excessive computing time, thereby solving the corresponding tasks more efficiently. However, quantum computing has only been proven to be effective in gene assembly tasks.

[0097] From the above description, it can be seen that the tools in the related art have the following deficiencies:

[0098] 1. The tools in the related art use two independent optimization processes to handle genotyping and gene assembly;

[0099] 2. Genotyping and gene assembly use heuristic algorithms, which makes it difficult to find the global optimal solution;

[0100] 3. The de novo typing and assembly of polyploids (i.e., without a reference genome) is difficult to handle, and the processing framework for polyploids is quite different from that for diploids;

[0101] 4. It can only perform typing and assembly on diploid low-depth sequencing data, but cannot process polyploid low-depth sequencing data.

[0102] Secondly, quantum computing in related technologies has the following deficiencies:

[0103] 1. It only targets gene assembly and cannot process genotyping at the same time;

[0104] 2. No other data information was used to assist in gene assembly.

[0105] Based on this, the embodiments of the present application provide a genotyping assembly method and device, an electronic device and a storage medium, so as to be able to perform genotyping and gene assembly simultaneously.

[0106] The genotyping assembly method provided in the examples of the present application is described below.

[0107] Reference Figure 1 In some embodiments, the genotyping assembly method provided in the embodiments of the present application includes but is not limited to steps S101 to S103.

[0108] Step S101, obtaining target comparison data from sequencing data of a target genome sequence, wherein the target comparison data is used to represent comparison data of any two target fragment sequences in the target genome sequence;

[0109] Step S102, constructing a path planning model according to the target alignment data and the target fragment sequence, wherein the path planning model is used to represent a path for typing and assembling the target fragment sequence according to the original sequence data;

[0110] Step S103, taking the original sequence data corresponding to the path minimization of the path planning model as the target sequence data, and performing typing and assembly on the target fragment sequence according to the target sequence data.

[0111] Steps S101 to S103 shown in the embodiment of the present application, a path planning model is constructed by target comparison data and target fragment sequence, and the path planning model can be used to describe the path of typing and assembling the target fragment sequence according to the original sequence data. It can be seen from this that solving the path planning model is actually a combinatorial optimization problem, and the corresponding target sequence data when the path planning model is minimized can be used as the optimal solution of the combinatorial optimization problem. Among them, the path planning model can be processed by path minimization based on the purpose of genotyping and the purpose of gene assembly. Therefore, when the target fragment sequence is typed and assembled according to the target sequence data, genotyping and gene assembly can be achieved simultaneously.

[0112] In step S101 of some embodiments, the target genome sequence refers to the genome sequence that needs to be genotyped and assembled. The target genome sequence is sequenced by a sequencing method such as nanopore protein to obtain sequencing data. The target genome sequence includes a target fragment sequence (i.e., target reads), and the sequencing data includes relevant information of the target fragment sequence. Based on an alignment tool (such as minimap2, minimap2 is a software tool for aligning DNA and RNA sequencing data), the sequencing data is analyzed and aligned to obtain target alignment data. Among them, the target alignment data is used to represent the alignment data of any two target fragment sequences.

[0113] Reference Figure 2 In some embodiments, any two target fragment sequences are taken as a target sequence pair, and step S101 includes but is not limited to steps S201 to S204.

[0114] Step S201, obtaining original alignment data from sequencing data, wherein the original alignment data is used to represent alignment data of any two original fragment sequences in the target genome sequence, the original alignment data includes an original alignment starting position and an original overlapping region length, and any two original fragment sequences are taken as an original sequence pair;

[0115] Step S202, taking any original fragment sequence in the original sequence pair as a comparison sequence, and if the original comparison start position is in the front region of the comparison sequence, taking the original sequence pair as a candidate sequence pair;

[0116] Step S203, screening the candidate sequence pairs according to the length of the original overlapping region to obtain the target sequence pairs;

[0117] Step S204: using the original alignment data corresponding to the target sequence pair as target alignment data.

[0118] It is understandable that multiple original fragment sequences can be obtained by sequence interruption of the target genome sequence. The target fragment sequence is a sequence obtained by screening the multiple original fragment sequences. The screening method is described below.

[0119] In step S201 of some embodiments, the original comparison data refers to the comparison data of any two original fragment sequences. The original comparison data can be obtained by analyzing and comparing the sequencing data based on the comparison tool. Taking the comparison tool as minimap2 as an example, the original comparison data can be obtained according to the following instructions:

[0120] minimap2-ck 28-for-only-Dt 8-x ava-pb input.fasta input.fasta>output.paf

[0121] Among them, "input" in input.fasta is the file name of the sequencing data, and the sequencing data can be in fasta format. "output" in output.paf is the file name of the original alignment data, and the original alignment data can be in paf format. It can be understood that, according to actual needs, the sequencing data can also be in fastq format, and the original alignment data can also be in bam format.

[0122] Reference Figure 3 In some embodiments, the original comparison data may be Figure 3 The table format is shown. The original alignment data includes the original alignment starting position (i.e., query_start and target_start) and the original overlapping area length (i.e., alignment_block_length). It can be understood that for any two original fragment sequences, one of the original fragment sequences is used as query reads and the other original fragment sequence is used as target reads, and the query reads and target reads can be combined into an original sequence pair. Among them, query_start is the original alignment starting position of the original fragment sequence queryreads, and target_start is the original alignment starting position of the original fragment sequence target reads. alignment_block_length can be regarded as the total number of bases in the alignment area between the original fragment sequence query reads and the original fragment sequence target reads.

[0123] In step S202 of some embodiments, the original fragment sequence query reads or the original fragment sequence target reads are used as the comparison sequence, the original sequence pairs whose original comparison start position is in the front region of the comparison sequence are retained, and the remaining original sequence pairs are filtered out. The retained original sequence pairs are used as candidate sequence pairs. It is understandable that the filtered original sequence pairs may be in the following situations:

[0124] In the first case, there is no overlapping region between the query reads and target reads in the original sequence pair;

[0125] In the second case, the original alignment start position is both in the front region of the original fragment sequence query reads and in the front region of the original fragment sequence target reads;

[0126] In the third case, the starting position of the original alignment is neither in the front region of the original fragment sequence query reads nor in the front region of the original fragment sequence target reads.

[0127] Among them, for the first case, the comparison tool can automatically filter it out.

[0128] Taking the third case as an example, refer to Figure 4 , the original sequence pair is X1-Y1, the position of point q1 is the beginning of the original fragment sequence Y1, L1 is the comparison area of ​​the original sequence pair X1-Y1, and the position of point q2 is the original comparison start position of the original fragment sequence X1. Figure 4 It can be seen that the original alignment starting position is neither in the front end region of the original fragment sequence X1 nor in the front end region of the original fragment sequence Y1.

[0129] Reference Figure 5 In some embodiments, the original alignment starting position includes a first starting position and a second starting position, and in step S202, "if the original alignment starting position is in the front end region of the alignment sequence, the original sequence pair is used as a candidate sequence pair" includes but is not limited to steps S501 to S502.

[0130] Step S501, comparing the first starting position and the second starting position with a preset position threshold respectively to obtain a position comparison result;

[0131] Step S502: If the position comparison result indicates that the first starting position or the second starting position is smaller than a preset position threshold, the original sequence pair is used as a candidate sequence pair.

[0132] In step S501 of some embodiments, query_start can be used as the first starting position and target_start can be used as the second starting position. Alternatively, target_start can be used as the first starting position and query_start can be used as the second starting position. The first starting position and the second starting position are respectively compared with the preset position threshold to obtain the corresponding position comparison result. It can be understood that the preset position threshold can be adaptively set according to actual needs, such as being set to 20, and this embodiment of the present application does not specifically limit this.

[0133] In step S502 of some embodiments, if the position comparison result indicates that the first starting position or the second starting position is less than the preset position threshold, it indicates that the first starting position or the second starting position is in the front region of the comparison sequence, and the corresponding original sequence pair is used as a candidate sequence pair.

[0134] It is understandable that, with the preset position thresholds of 20 and Figure 4For example, when the original alignment starting position corresponding to point q2 is greater than or equal to 20, it indicates that the original fragment sequence Y1 has greater than or equal to 20 bases between point q1 and point q2 that cannot be aligned with the original fragment sequence X1. At this time, this alignment is considered to be of low quality and needs to be filtered out.

[0135] In some embodiments, in step S203 to step S204, the candidate sequence pairs are further screened according to the length of the original overlapping region, and the candidate sequence pairs retained after the screening are used as the target sequence pairs. It can be understood that the original alignment data corresponding to the target sequence pairs are the target alignment data.

[0136] The advantage of screening the original sequence pairs by the original alignment starting position and the front region of the alignment sequence in steps S201 to S204 is that the original sequence pairs with low alignment quality and no help for gene assembly can be filtered out.

[0137] Reference Figure 6 In some embodiments, the original comparison data further includes a first sequence length and a second sequence length, and step S203 includes but is not limited to steps S601 to S604.

[0138] Step S601, calculating a multiple of a first sequence length according to a first preset multiple to obtain a first multiple length;

[0139] Step S602, calculating the multiple of the second sequence length according to the first preset multiple to obtain a second multiple length;

[0140] Step S603, comparing the original overlapping area length with the first multiple length and the second multiple length respectively to obtain a first length comparison result;

[0141] Step S604: if the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, a target sequence pair is obtained according to the candidate sequence pair.

[0142] It should be noted that, refer to Figure 3 , the original alignment data also includes the first sequence length (such as query_length) and the second sequence length (such as target_length).

[0143] In step S601 of some embodiments, the first preset multiple refers to pre-set data for calculating the multiple of the first sequence length, for example, it can be 0.3. It is understandable that the first preset multiple can be adaptively set according to actual needs, and this embodiment of the present application is not specifically limited to this. The first preset multiple is multiplied by the first sequence length (i.e., multiple calculation) to obtain the first multiple length.

[0144] In step S602 of some embodiments, similar to step S601, the first preset multiple is multiplied by the second sequence length to obtain the second multiple length.

[0145] In steps S603 to S604 of some embodiments, the original overlapping region length alignment_block_length corresponding to the candidate sequence pair is compared with the first multiple length and the second multiple length respectively, and if the first length comparison result obtained by comparison indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, the corresponding candidate sequence pair is retained and used as the target sequence pair. In other words, the candidate sequence pairs whose original region lengths are greater than both the first multiple length and the second multiple length are filtered out.

[0146] The benefit of step S601 to step S604 is that step S601 to step S604 is used to filter out candidate sequence pairs whose original overlapping region length is too short. For the filtered candidate sequence pairs, the two original fragment sequences included in the candidate sequence pair are most likely not adjacent sequences. Assuming that ABCDE... is the correct splicing order of a group of fragment sequences, the filtered candidate sequence pair may be AD, or AE, that is, an incorrect splicing order. Although such candidate sequence pairs will not affect the subsequent path minimization process, there are many such candidate sequence pairs. If these candidate sequence pairs are retained, it will affect the difficulty of path minimization processing of the path planning model. Therefore, in an embodiment of the present application, these candidate sequence pairs are filtered out to ensure that the target sequence data can be obtained according to the path minimization process in the future, so as to realize the typing and assembly of the target fragment sequence according to the target sequence data.

[0147] Reference Figure 7 In some embodiments, step S604 includes but is not limited to steps S701 to S705.

[0148] Step S701, if the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, the candidate sequence pair is used as the initial sequence pair;

[0149] Step S702, calculating the multiple of the first sequence length according to the second preset multiple to obtain a third multiple length, wherein the second preset multiple is greater than the first preset multiple;

[0150] Step S703, calculating the multiple of the second sequence length according to the second preset multiple to obtain a fourth multiple length;

[0151] Step S704, comparing the original overlapping area length with the third multiple length and the fourth multiple length respectively to obtain a second length comparison result;

[0152] Step S705: if the second length comparison result indicates that the length of the original overlapping region is greater than or equal to the third multiple length, and the length of the original overlapping region is greater than or equal to the fourth multiple length, the initial sequence pair is used as the target sequence pair.

[0153] In step S701 of some embodiments, further screening may be performed based on steps S601 to S604, and in this case, the candidate sequence pairs retained in step S604 are used as the initial sequence pairs.

[0154] In step S702 of some embodiments, the second preset multiple refers to pre-set data for calculating the multiple of the first sequence length, for example, a decimal greater than 0, less than 1, and close to 1 can be selected, such as selecting the second preset multiple as 0.9. It is understandable that the second preset multiple can be adaptively set according to actual needs, and this embodiment of the present application is not specifically limited to this. The second preset multiple is multiplied by the first sequence length (i.e., multiple calculation) to obtain a third multiple length.

[0155] In step S703 of some embodiments, similar to step S702, the second sequence length is multipled according to the second preset multiple to obtain a fourth multiple length;

[0156] In steps S704 to S705 of some embodiments, the original overlapping region length alignment_block_length corresponding to the initial sequence pair is compared with the third multiple length and the fourth multiple length respectively, and if the second length comparison result obtained by comparison indicates that the original overlapping region length is greater than or equal to the third multiple length, and the original overlapping region length is greater than or equal to the fourth multiple length, the corresponding initial sequence pair is retained and used as the target sequence pair. In other words, the initial sequence pairs whose original region lengths are less than both the third multiple length and the fourth multiple length are filtered out.

[0157] The benefit of steps S701 to S705 is that initial sequence pairs with original overlapping regions that are too long can be filtered out.

[0158] It is understandable that the embodiment of the present application does not specifically limit the order of filtering out sequence pairs whose original overlapping region lengths are too long and filtering out sequence pairs whose original overlapping region lengths are too short, that is, steps S701 to S705 can be performed first, and then steps S601 to S604 can be performed.

[0159] In step S102 of some embodiments, a path planning model for determining genotyping and gene assembly is constructed based on the target comparison data and the target fragment sequence. It is understandable that the path planning model is a model constructed based on the idea of ​​vehicle routing problem (VRP), and the path planning model is used to describe the path of splicing different target fragment sequences according to the original sequence data when performing genotyping and gene assembly. Different original sequence data correspond to different paths. Among them, VRP describes the problem of cargo delivery route planning in reality, specifically describing that there are K vehicles in the delivery fleet, and the goods in the warehouse are to be delivered to n customers located in different locations, and the distance between the receiving locations of each customer is known. The problem to be solved by VRP is how to assign n customers to K vehicles, and in what order each vehicle visits the assigned customers, so that the total driving distance of the K access paths corresponding to the K vehicles is the shortest. Correspondingly, the target fragment sequence can be used as a customer, the original sequence data can be used as the order in which the vehicle visits the corresponding customer, and the K paths corresponding to the K vehicles can be used as the K ploids (i.e., K haplotypes) that can be spliced.

[0160] Reference Figure 8 In some embodiments, step S102 includes but is not limited to steps S801 to S802.

[0161] Step S801, calculating and obtaining typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights;

[0162] Step S802, constructing a path planning model according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence number of the target fragment sequence, and the preset haplotype number.

[0163] In step S801 of some embodiments, the target sequence pairs obtained by screening according to the above steps and the target alignment data corresponding to the target sequence pairs can be constructed to obtain an alignment matrix. In the alignment matrix, the row index can be the number corresponding to the target fragment sequence at the front position in the target sequence pair, the column index can be the number corresponding to the target fragment sequence at the back position in the target sequence pair, and the matrix value can be the typing assembly prompt data D calculated according to the target alignment data and the preset weights. a,b. Specifically, if the original alignment starting position corresponding to a certain target fragment sequence is in the front region of the target fragment sequence, then the target fragment sequence is regarded as the target fragment sequence in the front position, that is, the number of the target fragment sequence can be used as the row index. Conversely, the target fragment sequence is regarded as the target fragment sequence in the back position, that is, the number of the target fragment sequence can be used as the column index. It can be seen that the method for determining whether the target fragment sequence is in the front position or the back position is similar to step S202, that is, if the original alignment starting position corresponding to the target fragment sequence is less than 20, it indicates that the target fragment sequence is in the front position, otherwise it is in the back position. For example, referring to Fig. 9 , for the target sequence pair X2-Y2, the length of the original alignment region corresponding to the target sequence pair is L2, point q3 is the original alignment start position corresponding to the target fragment sequence Y2, and point q4 is the original alignment start position corresponding to the target fragment sequence X2. Since point q3 is in the front region of the target fragment sequence Y2 and point q4 is in the non-front region of the target fragment sequence X2, the position of the target fragment sequence Y2 is in front and the position of the target fragment sequence X2 is in the back.

[0164] It is understandable that the numbers corresponding to the target fragment sequences can be obtained from Figure 3 The query_name and target_name in the query_name and target_name are obtained. It should be noted that Figure 3 The query_name and target_name in are only exemplary and do not represent the numbering method in the original alignment data. In actual situations, each original fragment sequence will be assigned a consecutive positive integer number starting from 1.

[0165] Classification assembly prompt data D a,b It is data set based on genotyping purpose and gene assembly purpose. Genotyping assembly prompts data D a,b It is used to prompt whether the target fragment sequence a and the target fragment sequence b need to be typed and assembled, that is, the target fragment sequence a and the target fragment sequence b can be guided to be typed and assembled. It is understandable that in the embodiment of the present application, the purpose of genotyping is to guide the target fragment sequence on the same homologous chromosome to be spliced, and the purpose of gene assembly is to guide the target fragment sequence from the same region of the chromosome or from the adjacent region of the chromosome to be spliced. The preset weight is a penalty parameter for the preset target comparison data. The specific value of the preset weight can be adaptively set according to actual needs, and the embodiment of the present application is not specifically limited to this.

[0166] Reference Fig.10 In some embodiments, the target comparison data includes the target overlapping region length and the target number of bases, the preset weights include a first weight and a second weight, and step S801 includes but is not limited to steps S1001 to S1004.

[0167] Step S1001, performing product calculation based on the first weight and the target base number to obtain first prompt data;

[0168] Step S1002, performing ratio calculation based on the target overlapping region length and the target base number to obtain base ratio data;

[0169] Step S1003, performing product calculation based on the second weight and the base ratio data to obtain second prompt data;

[0170] Step S1004, performing a summation calculation based on the second prompt data, the first prompt data and the target overlapping area length to obtain typing assembly prompt data.

[0171] It should be noted that the original alignment data corresponding to the target sequence pair is used as the target alignment data, so the original overlapping region length in the original alignment data can be used as the target overlapping region length. Figure 3 , the target alignment data also includes a target base number NM, which indicates the number of bases that cannot be aligned in the overlapping region. In some embodiments, the target alignment data also includes residue_matches, which indicates the number of bases that can be completely aligned in the overlapping region. The target overlapping region length alignment_block_length, the target base number NM and residue_matches have the following relationship:

[0172] alignment_block_length=NM+residue_matches

[0173] In step S1001 of some embodiments, the first weight is a preset penalty parameter, and the specific value of the first weight can be adaptively set according to actual needs, which is not specifically limited in the embodiments of the present application. For example, the first weight is set to p1, and the first weight is multiplied by the target base number to obtain the first prompt data p1*M a,b Among them, M a,b It represents the number of bases that cannot be aligned between the target fragment sequence a and the target fragment sequence b in the overlapping region (ie, NM of the target fragment sequence a and the target fragment sequence b).

[0174] In step S1002 to step S1003 of some embodiments, the ratio of the target overlapping region length to the target number of bases is calculated to obtain base ratio data M a,b / L a,b Among them, L a,bIndicates the length of the overlapping region between the target fragment sequence a and the target fragment sequence b (i.e., alignment_block_length). The second weight is a pre-set penalty parameter. The specific value of the second weight can be adaptively set according to actual needs, and this embodiment of the present application does not specifically limit this. For example, the second weight is set to p2, and the second weight is multiplied by the base ratio data to obtain the second prompt data p2*M a,b / L a,b .

[0175] In step S1004 of some embodiments, the second prompt data p2*M is calculated according to the following formula (1): a,b / L a,b , the first prompt data p1*M a,b The length of the overlapping area with the target is L a,b Perform sum calculation to obtain typing assembly prompt data D a,b .

[0176] D a,b =-L a,b +p1*M a,b +p2*M a,b / L a,b ......Formula (1)

[0177] Reference Fig.11 In some other embodiments, step S801 includes but is not limited to steps S1101 to S1104.

[0178] Step S1101, multiplying the target base number according to the preset weights to obtain base impact data;

[0179] Step S1102, performing difference calculation based on the target overlapping region length and the base influence data to obtain third prompt data;

[0180] Step S1103, performing a sum calculation based on the target overlapping region length and the base influence data to obtain fourth prompt data;

[0181] Step S1104, performing ratio calculation based on the third prompt data and the fourth prompt data to obtain typing and assembly prompt data.

[0182] In step S1101 of some embodiments, the preset weight p and the target number of bases M are a,b Perform product calculation to obtain base impact data p*M a,b .

[0183] In some embodiments, in step S1102 to step S1104, the target overlapping area length L a,b and the target base number M a,bPerform difference calculation to obtain the third prompt data L a,b -p*M a,b . Set the target overlapping area length L a,b and the target base number M a,b Perform sum calculation to obtain the fourth prompt data L a,b +p*M a,b According to the third prompt data L a,b -p*M a,b , the fourth prompt data L a,b +p*M a,b The typing assembly prompt data D is calculated as follows (2): a,b .

[0184]

[0185] It can be understood that whether it is the construction typing assembly prompt data D shown in formula (1) a,b The method is still the construction of typing assembly prompt data D shown in formula (2) a,b In addition to formula (1) and formula (2), other methods for constructing typing and assembly prompt data D can be selected based on the genotyping purpose and the gene assembly purpose. a,b From equations (1) and (2), we can see that when the target base number M a,b When it increases, it will affect the typing assembly prompt data D a,b Increase.

[0186] In step S802 of some embodiments, a path planning model as shown in the following formula (3) can be constructed based on the concept of VRP.

[0187]

[0188] Wherein, n represents the total number of customers in VRP, and represents the total number of target fragment sequences (i.e., the total number of sequences) in the embodiment of the present application. K represents the total number of vehicles in VRP, and the number of paths that need to be planned. In the embodiment of the present application, K represents the preset number of haplotypes, which means that K-ploid typing assembly is required (K≥1), which corresponds to the need to assemble K haplotypes. k represents the k-th vehicle in VRP, and represents the k-th haplotype in the embodiment of the present application. x i,a,k In VRP, it indicates whether the k-th vehicle visits customer a in step i (i is calculated from the departure of the warehouse, and the warehouse is step 0). In the embodiment of the present application, it indicates the typing assembly state of the target fragment sequence a, that is, whether the k-th haplotype needs to be spliced ​​with the target fragment sequence a in step i. If x i,a,k =1, otherwise x i,a,k = 0. Similarly, x i+1,b,kIn the embodiment of the present application, it indicates whether the k-th vehicle visits customer b in the i+1 step. i,a,k =1 and x i+1,b,k =1, indicating that the original sequence data is a→b. a,b In VRP, it represents the distance between customer a and customer b.

[0189] From the above description, we can see that the path planning model is a combinatorial optimization problem about the splicing order of the target fragment sequence. The solution of the path planning model corresponds to different original sequence data, and the target sequence data is the optimal solution of the path planning model. The optimal solution of the path planning model is obtained when the value of formula (3) is the smallest (that is, the path of the path planning model is minimized), that is, the typing assembly prompt data D a,b The smaller the value of is, the more conducive it is to splice the target fragment sequence a with the target fragment sequence b. Therefore, formula (1) and formula (2) construct the typing assembly prompt data D a,b The method can be used to determine the target base number M a,b A certain penalty is imposed on the target sequence with a high overlap area length L a,b Longer, higher quality alignment (such as target base number M a,b low) target sequence pairs, even if the target base number M a,b The high target sequence does not correspond to the optimal solution of the path planning model, and the target overlap area length L a,b The optimal solution of the target sequence with long length and high alignment quality is close to the path planning model. In this way, the genotyping purpose and the gene assembly purpose can be achieved at the same time. Specifically, there is a difference in SNP sites between two target fragment sequences from different homologous chromosomes, and the target base numbers M corresponding to these two target fragment sequences are a,b The value is high, based on the typing assembly prompt data D a,b The construction method and the method for solving the optimal solution of the path planning model will punish the two target fragment sequences for not splicing. In other words, the embodiment of the present application rewards and guides the two target fragment sequences from the same homologous chromosome to be spliced ​​to achieve the purpose of genotyping.

[0190] Secondly, for two target fragment sequences from the same homologous chromosome, if the two target fragment sequences come from different regions of the chromosome, the target base numbers M corresponding to the two target fragment sequences are a,b It will also be higher than the target base number M corresponding to other situations a,b. Among them, other situations refer to the same region or adjacent regions of the source chromosome. Similarly, two target fragment sequences from different regions of the chromosome will be punished for not being spliced. In other words, the embodiment of the present application rewards and guides the splicing of two target fragment sequences from the same homologous chromosome and from the same region or adjacent regions of the chromosome to achieve the purpose of gene assembly.

[0191] It is understandable that in the comparison matrix, there may be matrix values ​​(i.e., D a,b ) is less than 0 or equal to 0. In this case, in order to solve the path planning model, the minimum value min(D) of all matrix values ​​is determined. In the comparison matrix, except for the matrix values ​​on the diagonal of the comparison matrix, all other matrix values ​​are calculated with the minimum value min(D) to ensure that the minimum matrix value in the comparison matrix after the difference calculation is 0, thereby meeting the solution requirements of the path planning model.

[0192] The benefit of step S801 to step S802 is that the path planning model constructed based on the VRP concept can perform genotyping and gene assembly at the same time.

[0193] In step S103 of some embodiments, VRP is a combinatorial optimization problem, and the goal of VRP is to minimize the total driving distance of K access paths corresponding to K vehicles. Correspondingly, the goal of the path planning model is to find the original sequence data that minimizes the path. The original sequence data corresponding to the path minimization of the path planning model is used as the target sequence data, and the corresponding target fragment sequence is spliced ​​according to the target sequence data to obtain the final DNA sequence.

[0194] Reference Fig.12 In some embodiments, minimization processing of the path planning model includes but is not limited to steps S1201 to S1202.

[0195] Step S1201, constructing constraint conditions, wherein the constraint conditions include that the number of typing and assembly for each target fragment sequence is equal to 1, only one target fragment sequence can be typed and assembled at a time, and the typing and assembly starts from a preset start and end sequence;

[0196] Step S1202: performing path minimization processing on the path planning model according to the constraint conditions.

[0197] In step S1201 of some embodiments, when solving the VRP problem, it is necessary to combine the constraint conditions shown in equation (6) of equation (4) below.

[0198]

[0199]

[0200]

[0201] Among them, in the VRP problem, the constraint corresponding to formula (4) means that each vehicle must start from the warehouse, the constraint corresponding to formula (5) means that each vehicle only visits one customer at a time, and the constraint corresponding to formula (6) means that each customer is visited and only visited once. In the embodiment of the present application, the above constraints are transformed into starting from a preset start and end sequence (equivalent to the warehouse in the VPR problem) to perform typing assembly, only one target fragment sequence can be typed and assembled at a time, and the number of times each target fragment sequence is typed and assembled is equal to 1.

[0202] Reference Fig.13 In some embodiments, the method for constructing a preset start and end sequence includes but is not limited to steps S1301 to S1302.

[0203] Step S1301, constructing an original start and end sequence;

[0204] Step S1302 : compare the original start-end sequence with the target fragment sequence. If the original start-end sequence and the target fragment sequence do not include an overlapping region, the original start-end sequence is used as a preset start-end sequence.

[0205] It should be noted that in the VRP problem, the vehicle must start from the warehouse and return to the warehouse to form a circular loop. However, there is no circular DNA in gene assembly, so a virtual sequence fragment (i.e., preset start and end sequences) needs to be set to form a virtual loop. When determining the final DNA, it can be disconnected at the preset start and end sequences, that is, only the splicing order of the real target fragment sequence is concerned.

[0206] In step S1301 of some embodiments, an original start and end sequence is constructed, and the numbers of the original start and end sequences are set to 0.

[0207] In step S1302 of some embodiments, the original start-end sequence is compared with any target fragment sequence. If the original start-end sequence has no overlapping area with any target fragment sequence and the NM value is equal to 0, the original start-end sequence is used as the preset start-end sequence. In other words, the matrix value of the preset start-end sequence in the comparison matrix should satisfy the condition shown in formula (7).

[0208]

[0209] Reference Fig.14 In some embodiments, the constraint condition also includes a chromosome homology state indicating that any two target fragment sequences originate from the same homologous chromosome, and step S802 includes but is not limited to steps S1401 to S1403.

[0210] Step S1401, constructing an initial model according to typing assembly prompt data, typing assembly status of the target fragment sequence, total sequence quantity of the target fragment sequence, and a preset number of haplotypes;

[0211] Step S1402, determining the chromosome homology status of any two target fragment sequences;

[0212] Step S1403, constructing a path planning model according to the chromosome homology state and the initial model;

[0213] It should be noted that the typing assembly can also be assisted based on Hi-C (High-throughput chromosome conformation capture, a technology derived from chromatin conformation capture combined with high-throughput sequencing, which can realize the capture of interactions between chromosome fragments within the whole genome) sequencing data. Therefore, in step S1401, the model (such as the model shown in formula (3)) constructed according to the typing assembly prompt data, the typing assembly state of the target fragment sequence, the total sequence amount of the target fragment sequence, and the preset number of haplotypes is used as an intermediate model, that is, as an initial model. In steps S1402 and S1403, the initial model is combined with the Hi-C sequencing data to obtain the final path planning model.

[0214] In step S1402 of some embodiments, the Hi-C sequencing data obtains a pair of reads, and it can be determined that the pair of reads is most likely from the same homologous chromosome. It can be seen that the Hi-C sequencing data can be used to assist genotyping. The chromosome homology state of any two target fragment sequences (i.e., target sequence pairs) is determined based on the Hi-C sequencing data. The chromosome homology state is used to indicate that the corresponding target sequence pairs are from the same homologous chromosome, or from different homologous chromosomes.

[0215] In step S1403 of some embodiments, the chromosome homology state is combined with the initial model to construct a final path planning model as shown in the following formula (8).

[0216]

[0217] Among them, Q VRP It is the binary quadratic unconstrained optimization form of the initial model, a and b are chromosome homology states representing two target fragment sequences from the same homologous chromosome. PD represents a collection of multiple target sequence pairs, and the chromosome homology states of the target sequence pairs in PD all represent that they are from the same homologous chromosome.

[0218] Among them, Q VRP The specific form is shown as follows:

[0219]

[0220] It is understandable that when the final path planning model is determined in combination with Hi-C sequencing data, the constraints for solving the path planning model should also be increased accordingly. For example, in addition to the constraints corresponding to equations (4) to (6), the following two constraints need to be added (using the VRP problem for illustration):

[0221] 1. In terms of vehicle selection, the pickup vehicle and the delivery vehicle must be the same vehicle;

[0222] 2. In terms of the order of customer visits, the vehicle must first visit the customer who wants to pick up the goods, and then visit the customer who wants to deliver the goods.

[0223] In the typing and assembly problem of the embodiment of the present application, usually only the first constraint can be achieved, because according to the Hi-C sequencing data, it can only be determined that the two target fragment sequences are derived from the same homologous chromosome, but it is not possible to obtain information about which target fragment sequence is spliced ​​first and which target fragment sequence is spliced ​​later. Therefore, the newly added constraint is that the chromosome homology state indicates that any two target fragment sequences (i.e., target sequence pairs) are derived from the same homologous chromosome.

[0224] The benefit of step S1401 to step S1403 is that Hi-C sequencing data can be combined to assist typing and assembly, thereby improving the accuracy of typing and assembly.

[0225] In step S1202 of some embodiments, the optimal solution of the path planning model is obtained according to the constraint conditions, that is, the path minimization processing is performed on the path planning model to obtain the target sequence data corresponding to the K haplotypes. It can be understood that in the VRP problem, the optimal solution can be which customers a certain vehicle needs to visit in what order. Corresponding to the embodiment of the present application, the optimal solution of the path planning model can be to splice the corresponding target fragment sequence according to the target sequence data to obtain one of the K haplotypes.

[0226] In some embodiments, the path planning model can be solved based on a quantum algorithm. When using a quantum algorithm, the problem to be solved needs to be encoded on a quantum annealing chip (such as a D-Wave quantum annealing chip). Since the form shown in equation (3) cannot be directly encoded, Q is used when performing quantum computing. VRP It is understandable that, according to actual needs, other quantum computing hardware based on quantum annealing technology, classical optimization solvers (such as OR-Tools), coherent Ising machines (CIM), simulated quantum inspired bifurcation machines (SQBM), programmable quantum computing frameworks, etc. can also be used to solve the path planning model.

[0227] Reference Fig.15 In some embodiments, step S1202 includes but is not limited to steps S1501 to S1502.

[0228] Step S1501, obtaining annealing time and annealing times;

[0229] Step S1502, performing path minimization processing on the path planning model according to the preset quantum annealing model, annealing time, annealing number, and constraint conditions.

[0230] In step S1501 of some embodiments, the calculation parameters of the quantum annealing chip, namely the annealing duration and the number of annealing times, are obtained. The annealing duration can be set in advance according to the path planning model. The annealing duration is not the longer the better. Excessive annealing duration will cause the quantum system inside the quantum annealing chip to decoherent, resulting in a loss of credibility in the final result. Based on the performance of the quantum annealing chip, the annealing duration can be set to any value between 20μs and 40μs, and the number of annealing times can be set to 1000 times. The annealing duration and the number of annealing times can also be adaptively set according to actual needs, which is not specifically limited in the embodiments of the present application.

[0231] In step S1502 of some embodiments, the quantum annealing chip (i.e., quantum annealing model) is set according to the annealing duration and the number of annealing times, and the interface of the quantum annealing chip is called to pass the path planning model into the quantum annealing chip as a dictionary type, so as to solve the path planning model based on the quantum annealing chip, that is, perform path minimization processing on the path planning model.

[0232] It can be understood that, according to the target sequence data obtained by path minimization processing, the corresponding target fragment sequences are spliced, and K monomeric DNA sequences can be spliced ​​according to multiple target sequence data.

[0233] In a specific embodiment, referring to Fig.16 , take the simulated triploid short DNA sequence as an example (that is, the triploid short DNA sequence is used as the target genome sequence in the above embodiment). The three short DNA sequences (including homologous chromosome 1 DNA, homologous chromosome 2 DNA, and homologous chromosome 3 DNA) are broken into fragments as model reads (that is, the original fragment sequence). Assume that the first weight p1=5, the second weight p2=0, the species ploidy K=3, the annealing time is 30μs, and the number of annealing times is 1000 times. According to the description of the above embodiment, three target sequence data can be obtained. According to the three target sequence data, the corresponding target fragment sequences are typed and assembled to obtain the following Fig.17The three DNA sequences shown (including the DNA sequence corresponding to the assembly result of haplotype 1, the DNA sequence corresponding to the assembly result of haplotype 2, and the DNA sequence corresponding to the assembly result of haplotype 3). Fig.16 and Fig.17 It can be seen from the comparison that the assembly result of the embodiment of the present application is correct, which proves the feasibility of the embodiment of the present application. Fig.18 As shown, in the above example, a single annealing takes 30 μs. Fig.18 The 193.000ms in the figure is the total annealing time of the quantum annealing chip. The classical combinatorial optimization solver Gurobi is used to solve this example, which takes 0.3s. By comparing the time consumption, it can be seen that the quantum annealing chip provides three orders of magnitude of computational acceleration, demonstrating the efficiency of quantum computing in small-scale analytical assembly problems.

[0234] The genotyping assembly method provided in the embodiment of the present application realizes the grouping assembly of any ploid, combines genotyping and gene assembly into one, and simplifies the gene analysis process. By setting a suitable preset weight (i.e., penalty parameter), annealing time and annealing number, through the solution space exploration ability of a quantum computer, the assembly quality of the heuristic bioinformatics tool in the relevant technology is obtained, the typing accuracy and assembly continuity are improved, and the DNA sequence of the typing assembly is made more accurate. In addition, based on the manipulation time of quantum computing nanoseconds, the embodiment of the present application shortens the calculation time of polyploid typing assembly. In addition, the typing assembly process in the relevant technology is run on a large-scale computing cluster, and the maintenance and energy consumption control of the cluster are crucial. This fundamental reason is that the operation of a classical computer is an irreversible calculation, and irreversible calculation can bring about a thermal effect. The embodiment of the present application uses the characteristics of reversible calculations run by a quantum computer to significantly reduce the heat and energy consumption generated when performing a set of polyploid typing assembly processes.

[0235] The genotyping assembly device provided in the embodiments of the present application is described below.

[0236] Reference Fig.19 The present application also provides a genotyping assembly device, which includes:

[0237] The comparison data acquisition module 1901 is used to obtain target comparison data from the sequencing data of the target genome sequence, wherein the target comparison data is used to represent the comparison data of any two target fragment sequences in the target genome sequence;

[0238] A path planning model building module 1902 is used to build a path planning model based on the target alignment data and the target fragment sequence, wherein the path planning model is used to represent the path for typing and assembling the target fragment sequence based on the original sequence data;

[0239] The typing and assembling module 1903 is used to use the original sequence data corresponding to the path minimization of the path planning model as the target sequence data, and to perform typing and assembly on the target fragment sequence according to the target sequence data.

[0240] It can be seen that the contents of the above-mentioned genotyping assembly method embodiment are all applicable to the embodiment of the present genotyping assembly device, the functions specifically implemented by the present genotyping assembly device embodiment are the same as those of the above-mentioned genotyping assembly method embodiment, and the beneficial effects achieved are also the same as those achieved by the above-mentioned genotyping assembly method embodiment.

[0241] Reference Fig. 20 , Fig. 20 The hardware structure of an electronic device of another embodiment is illustrated, and the electronic device includes:

[0242] The processor 2001 may be implemented by a general-purpose CPU (Central Processing Unit), a microprocessor, an application-specific integrated circuit (Application Specific Integrated Circuit, ASIC), or one or more integrated circuits, and is used to execute relevant programs to implement the technical solutions provided in the embodiments of the present application;

[0243] The memory 2002 can be implemented in the form of a read-only memory (ROM), a static storage device, a dynamic storage device, or a random access memory (RAM). The memory 2002 can store an operating system and other application programs. When the technical solution provided in the embodiments of this specification is implemented by software or firmware, the relevant program code is stored in the memory 2002, and the processor 2001 is used to call and execute the genotyping assembly method of the embodiment of this application;

[0244] Input / output interface 2003, used to implement information input and output;

[0245] Communication interface 2004, used to realize communication interaction between the device and other devices, which can be realized through wired mode (such as USB, network cable, etc.) or wireless mode (such as mobile network, WIFI, Bluetooth, etc.);

[0246] Bus 2005 , which transmits information between various components of the device (e.g., processor 2001 , memory 2002 , input / output interface 2003 , and communication interface 2004 );

[0247] The processor 2001 , the memory 2002 , the input / output interface 2003 and the communication interface 2004 are connected to each other in communication within the device via the bus 2005 .

[0248] The embodiment of the present application also provides a computer program product, which includes a computer program. The processor of the computer device reads and executes the computer program, so that the computer device executes the above-mentioned genotyping assembly method.

[0249] The terms "first", "second", "third", "fourth", etc. (if any) in the specification of the present disclosure and the above-mentioned drawings are used to distinguish similar objects, and are not necessarily used to describe a specific order or sequence. It should be understood that the data used in this way can be interchangeable where appropriate, so that the embodiments of the present disclosure described herein can, for example, be implemented in an order other than those illustrated or described herein. In addition, the terms "comprises" and "comprising" and any variations thereof are intended to cover non-exclusive inclusions, for example, a process, method, system, product or device that includes a series of steps or units is not necessarily limited to those steps or units that are clearly listed, but may include other steps or units that are not clearly listed or inherent to these processes, methods, products or devices.

[0250] It should be understood that in the present disclosure, "at least one (item)" means one or more, and "plurality" means two or more. "And / or" is used to describe the association relationship of associated objects, indicating that three relationships may exist. For example, "A and / or B" can mean: only A exists, only B exists, and A and B exist at the same time, where A and B can be singular or plural. The character " / " generally indicates that the objects associated before and after are in an "or" relationship. "At least one of the following" or similar expressions refers to any combination of these items, including any combination of single or plural items. For example, at least one of a, b or c can mean: a, b, c, "a and b", "a and c", "b and c", or "a and b and c", where a, b, c can be single or multiple.

[0251] It should be understood that in the description of the embodiments of the present application, the meaning of multiple (or multiple items) is more than two, greater than, less than, exceed, etc. are understood to not include the number, and above, below, within, etc. are understood to include the number.

[0252] In the several embodiments provided in the present disclosure, it should be understood that the disclosed systems, devices and methods can be implemented in other ways. For example, the device embodiments described above are only schematic. For example, the division of units is only a logical function division. There may be other division methods in actual implementation. For example, multiple units or components can be combined or integrated into another system, or some features can be ignored or not executed. Another point is that the mutual coupling or direct coupling or communication connection shown or discussed can be an indirect coupling or communication connection through some interfaces, devices or units, which can be electrical, mechanical or other forms.

[0253] The units described as separate components may or may not be physically separated, and the components shown as units may or may not be physical units, that is, they may be located in one place or distributed on multiple network units. Some or all of the units may be selected according to actual needs to achieve the purpose of the solution of this embodiment.

[0254] In addition, each functional unit in each embodiment of the present disclosure may be integrated into one processing unit, or each unit may exist physically separately, or two or more units may be integrated into one unit. The above-mentioned integrated unit may be implemented in the form of hardware or in the form of software functional units.

[0255] If the integrated unit is implemented in the form of a software functional unit and sold or used as an independent product, it can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present disclosure is essentially or the part that contributes to the prior art or all or part of the technical solution can be embodied in the form of a software product, and the computer software product is stored in a storage medium, including a number of instructions for a computer device (which can be a personal computer, a server, or a network device, etc.) to perform all or part of the steps of the various embodiments of the present disclosure. The aforementioned storage medium includes: U disk, mobile hard disk, read-only memory (Read-Only Memory, referred to as ROM), random access memory (Random Access Memory, referred to as RAM), disk or optical disk and other media that can store program codes.

[0256] It should also be understood that the various implementation methods provided in the embodiments of the present application can be combined arbitrarily to achieve different technical effects.

[0257] The above is a specific description of the implementation methods of the present disclosure, but the present disclosure is not limited to the above implementation methods. Technical personnel familiar with the art can also make various equivalent modifications or substitutions without violating the spirit of the present disclosure. These equivalent modifications or substitutions are all included in the scope defined by the claims of the present disclosure.

Claims

1. A genotyping assembly method, characterized in that: The method comprises: Obtaining target comparison data from the sequencing data of the target genome sequence, wherein the target comparison data is used to represent the comparison data of any two target fragment sequences in the target genome sequence; Constructing a path planning model according to the target comparison data and the target fragment sequence, wherein the path planning model is used to represent a path for typing and assembling the target fragment sequence according to the original sequence data; The original sequence data corresponding to the path minimization of the path planning model is used as the target sequence data, and the target fragment sequence is typed and assembled according to the target sequence data.

2. The method according to claim 1, characterized in that Taking any two of the target fragment sequences as a target sequence pair; The step of obtaining target comparison data from sequencing data of a target genome sequence comprises: Obtaining original comparison data from the sequencing data, wherein the original comparison data is used to represent the comparison data of any two original fragment sequences in the target genome sequence, the original comparison data includes the original comparison starting position and the original overlapping region length, and taking any two original fragment sequences as an original sequence pair; Taking any original fragment sequence in the original sequence pair as a comparison sequence, if the original comparison start position is in the front region of the comparison sequence, taking the original sequence pair as a candidate sequence pair; Screening the candidate sequence pairs according to the length of the original overlapping region to obtain the target sequence pairs; The original alignment data corresponding to the target sequence pair is used as the target alignment data.

3. The method according to claim 2, characterized in that The original comparison starting position includes a first starting position and a second starting position, and if the original comparison starting position is in the front end region of the comparison sequence, taking the original sequence pair as a candidate sequence pair includes: Compare the first starting position and the second starting position with a preset position threshold value respectively to obtain a position comparison result; If the position comparison result indicates that the first starting position or the second starting position is smaller than the preset position threshold, the original sequence pair is used as a candidate sequence pair.

4. The method according to claim 2, characterized in that The original comparison data also includes a first sequence length and a second sequence length; The step of screening the candidate sequence pairs according to the length of the original overlapping region to obtain the target sequence pairs comprises: Calculate the first sequence length by multiples according to a first preset multiple to obtain a first multiple length; Calculate the multiple of the second sequence length according to the first preset multiple to obtain a second multiple length; Compare the original overlapping area length with the first multiple length and the second multiple length respectively to obtain a first length comparison result; If the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, the target sequence pair is obtained according to the candidate sequence pair.

5. The method according to claim 4, characterized in that If the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, obtaining the target sequence pair according to the candidate sequence pair includes: If the first length comparison result indicates that the original overlapping region length is less than or equal to the first multiple length, and the original overlapping region length is less than or equal to the second multiple length, taking the candidate sequence pair as an initial sequence pair; Calculating the first sequence length by multiples according to a second preset multiple to obtain a third multiple length, wherein the second preset multiple is greater than the first preset multiple; Calculate the multiple of the second sequence length according to the second preset multiple to obtain a fourth multiple length; Compare the original overlapping area length with the third multiple length and the fourth multiple length respectively to obtain a second length comparison result; If the second length comparison result indicates that the original overlapping region length is greater than or equal to the third multiple length, and the original overlapping region length is greater than or equal to the fourth multiple length, the initial sequence pair is used as the target sequence pair.

6. The method according to claim 1, characterized in that The step of constructing a path planning model according to the target comparison data and the target fragment sequence includes: Calculate the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights; The path planning model is constructed according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence quantity of the target fragment sequence, and the preset haplotype quantity.

7. The method according to claim 6, characterized in that The target comparison data includes the target overlapping region length and the target base number, and the preset weights include a first weight and a second weight; The step of calculating the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights includes: Perform product calculation based on the first weight and the target base number to obtain first prompt data; Calculating the ratio of the target overlapping region length to the target base number to obtain base ratio data; Perform product calculation based on the second weight and the base ratio data to obtain second prompt data; The typing assembly prompt data is obtained by performing a sum calculation based on the second prompt data, the first prompt data and the target overlapping area length.

8. The method according to claim 6, characterized in that The target comparison data includes the target overlapping region length and the target base number; The step of calculating the typing and assembly prompt data of any two target fragment sequences according to the target comparison data and preset weights includes: Perform product calculation on the target number of bases according to the preset weights to obtain base impact data; Performing difference calculation based on the target overlapping region length and the base impact data to obtain third prompt data; A fourth prompt data is obtained by performing a sum calculation based on the target overlapping region length and the base influence data; The typing and assembly prompt data are obtained by performing a ratio calculation based on the third prompt data and the fourth prompt data.

9. The method according to claim 6, characterized in that Before the original sequence data corresponding to the path minimization of the path planning model is used as the target sequence data, the method further includes performing path minimization processing on the path planning model, specifically including: Constructing constraint conditions, wherein the constraint conditions include that the number of typing and assembling of each target fragment sequence is equal to 1, only one target fragment sequence can be typed and assembled at a time, and the typing and assembly starts from a preset start and end sequence; The path planning model is subjected to path minimization processing according to the constraint conditions.

10. The method according to claim 9, characterized in that Before constructing the constraint data, the method further includes constructing a preset start and end sequence, specifically including: Construct the original start-end sequence; The original start-end sequence is compared with the target fragment sequence. If the original start-end sequence and the target fragment sequence do not include an overlapping region, the original start-end sequence is used as the preset start-end sequence.

11. The method according to claim 9, characterized in that The step of constructing the path planning model according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence number of the target fragment sequence, and the preset number of haplotypes includes: Constructing an initial model according to the typing assembly prompt data, the typing assembly status of the target fragment sequence, the total sequence quantity of the target fragment sequence, and the preset number of haplotypes; Determining the chromosome homology status of any two target fragment sequences; Constructing the path planning model according to the homologous state of the chromosome and the initial model; The constraint condition also includes that the chromosome homology state indicates that any two target fragment sequences originate from the same homologous chromosome.

12. The method according to claim 9, characterized in that The performing path minimization processing on the path planning model according to the constraint condition comprises: Get the annealing time and annealing times; The path planning model is subjected to path minimization processing according to a preset quantum annealing model, the annealing duration, the number of annealing times, and the constraint conditions.

13. A genotyping assembly device, characterized in that: The device comprises: A comparison data acquisition module, used to obtain target comparison data from the sequencing data of the target genome sequence, wherein the target comparison data is used to represent the comparison data of any two target fragment sequences in the target genome sequence; A path planning model construction module, used to construct a path planning model according to the target comparison data and the target fragment sequence, wherein the path planning model is used to represent a path for typing and assembling the target fragment sequence according to the original sequence data; The typing and assembly module is used to use the original sequence data corresponding to the path minimization of the path planning model as the target sequence data, and to perform typing and assembly on the target fragment sequence according to the target sequence data.

14. An electronic device comprising a memory and a processor, wherein the memory stores a computer program, wherein: When the processor executes the computer program, the genotyping assembly method according to any one of claims 1 to 12 is implemented.

15. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by a processor, the genotyping assembly method according to any one of claims 1 to 12 is implemented.

16. A computer program product, comprising a computer program, wherein the computer program is read and executed by a processor of a computer device, so that the computer device executes the genotyping assembly method according to any one of claims 1 to 12.

Images