A graph-based multiple sequence alignment method and system

Abstract

The present invention discloses a graph-based multiple sequence alignment method, system, device and computer-readable storage medium, the method comprising: obtaining a set of pan-genome sequence data, mapping the pan-genome sequence data, and obtaining a pan-genome colored graph; marking and obtaining the characteristics of the access status of a single node of the colored graph, traversing the colored graph to obtain a cSupB data model after decomposing the colored graph; obtaining the cSupB data model after the colored graph is decomposed, extracting the offset value data characteristics of the nodes in the cSupB data model; performing a first preprocessing on the cSupB data model based on the offset value data characteristics, and obtaining the initial alignment result of the cSupB data model. The present invention overcomes the problem that the traditional multiple sequence alignment method cannot meet the current situation when facing a large number of sequenced sample genomes.

Description

Technical Field

[0001] The present invention belongs to the field of medical technology, and in particular relates to a graph-based multiple sequence alignment method and system thereof. Background Art

[0002] The development of life sciences, medicine and other fields is closely related to the application of sequencing technology. However, due to sequencing technology, sequencing costs and even computing costs, many genome studies have many problems, such as over-reliance on reference genomes. At present, reference genomes occupy a very important position in many fields. In almost all studies involving genomes, the first thing people do is to construct a reference genome for the species being studied, and then conduct different follow-up studies based on the reference genome, such as comparing the data of other newly sequenced individuals of the species with the reference genome to find differences. This method is the basis for seeking the origin of disease genes in human genomics. However, the biggest disadvantage of the reference genome-based method is the omission problem, because a single genome obviously cannot contain all the information of the genome. In the context of today's large number of species and individuals that can be widely sequenced, taking humans as an example, if the traditional reference genome research method is still used, at least 10% of the human genome sequence information will be omitted in the reference genome.

[0003] In recent years, with the development of sequencing technology, the quality of individual genome assembly has been continuously improved, and the reduction of sequencing costs has also led to an increasing number of sequencing. Taking humans as an example, the assembly quality of sequenced sample genomes is already comparable to GRCh38. There are many available genome assembly results. I believe that this number will continue to increase in the future, not only for humans, but also for other species. We are gradually entering the population genome era from the genome era. The advent of the population genome era has brought a large amount of genomic data and unprecedented research opportunities, but also put forward new requirements and challenges for bioinformatics analysis methods. For example, how to effectively organize large-scale population genome data and conduct subsequent analysis (such as phylogenetic analysis) is a problem that researchers urgently need to solve.

[0004] At present, the basic process commonly used for phylogenetic analysis of a group of nucleic acid sequences includes multiple sequence alignment, conserved region selection, saturation analysis, evolutionary model selection, and developmental tree construction. However, with the increase in the number of sequences, it is obviously no longer feasible to directly align all sequences. Most of the traditional commonly used multiple sequence alignment algorithms rely on progressive alignment algorithms. This algorithm aligns the input sequences one by one according to the distance of the relationship according to the pre-calculated guide tree, and then merges them into the final model. However, due to the high similarity of population genomes, long sequences, and large sample sizes, the estimation of guide trees and progressive alignment will inevitably lead to many redundant operations, greatly increasing the complexity of the calculation, and the cost of the entire alignment process is high. In addition, if the method of pre-processing group alignment is adopted, it will also consume a lot of time and computing resources, and directly ignore the relationship between individual sequences between groups.

[0005] In order to solve the problem that the traditional multiple sequence alignment method cannot meet the current situation when facing a large number of sequenced sample genomes, the present invention provides a graph-based multiple sequence alignment method and system. Summary of the invention

[0006] In order to overcome the problems mentioned in the background technology, the present invention provides a graph-based multiple sequence alignment method and system thereof.

[0007] A graph-based multiple sequence alignment method, comprising:

[0008] Acquire a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome;

[0009] Mark and obtain the characteristics of the access status of a single node of the colored graph, traverse the colored graph to obtain a cSupB data model after decomposing the colored graph;

[0010] Obtaining a cSupB data model after the colored graph is decomposed, and extracting offset value data features of nodes in the cSupB data model;

[0011] Performing a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model;

[0012] Optionally, the characteristics of the access state of a single node of the colored graph include an unvisited state, a semi-visited state, an accessible state and a visited state;

[0013] Optionally, the unvisited state is that no entry point has been visited; the semi-visited state is that at least one entry point has been visited and at least one entry point has not been visited; the accessible state is that all entry points have been visited and the node itself is in a state that can be accessed; the visited state is that all entry points of the node have been visited and the node itself has also been visited;

[0014] Optionally, the method of traversing the colored graph includes but is not limited to a post-order traversal method;

[0015] Optionally, the cSupB data model is: in a colored graph G = (V, E, C), V(G), V(E) and V(C) are the point set, edge set and color set of the graph G respectively. For any color set G1=(V1,E1,C1) is a subgraph of graph G, satisfying that for any node u i ∈V1, For two different points s and t, it is called<s,t,C1> A colored SuperBubble; s is called the source node, and t is the sink node;

[0016] Optionally, the colored graph is a directed, acyclic, non-degenerate base graph with a unique starting point and end point.

[0017] The first preprocessing includes: defining the cSupB data model as an s×l cSupB alignment matrix based on the offset value data feature, traversing the cSupB alignment matrix, adding the base information of each node to the corresponding position of the cSupB alignment matrix, and obtaining an initial alignment result of the cSupB data model;

[0018] Optionally, in the cSupB comparison matrix, s is the number of samples of the target cSupB data model, l=pos(source)-pos(sink)+1, and l is the maximum length of the target cSupB data model;

[0019] Optionally, the traversal process includes but is not limited to the following traversal methods: a post-order traversal method;

[0020] Optionally, the base information of each node is the last base of the k-mer.

[0021] Screening the cSupB data model to identify a target cSupB data model;

[0022] Optionally, the screening process includes: extracting a root cSupB data model of level=1 in the cSupB data model; the root cSupB data model is the target cSupB data model;

[0023] Optionally, the target cSupB data model is preprocessed for the first time to obtain an initial comparison result.

[0024] The multiple sequence alignment method further comprises: analyzing the initial alignment result to obtain a data analysis result of the initial alignment result, comparing the relationship between the data analysis result and the adjustment standard, determining whether to perform adjustment processing on the initial alignment result according to the relationship; performing or not performing adjustment processing, and outputting a final alignment result;

[0025] Optionally, the process of analyzing the initial alignment result includes: identifying the variation type of each aligned site in the initial alignment result, analyzing the initial alignment result column by column to obtain a hash table pos2vartype, traversing the hash table pos2vartype to obtain a data analysis result;

[0026] Optionally, the variation types include: no gap, 1 base, no variation; no gap, 2 bases, SNP; no gap, 3 or 4 bases, SNP; with gap, 1 base, Indel; with gap, more than 1 base, Indel+SNP;

[0027] Optionally, the key of the hash table pos2vartype is the matrix position, and the value is the variation type;

[0028] Optionally, the method for traversing the hash table pos2vartype is: traversing in sequence according to the size of the matrix position, encountering a gap area during the traversal, the output length is the gap length, calculating the misalignment score of the m bp sequences on both sides of the gap area, and the misalignment score is the data analysis result;

[0029] Optionally, the variation type of the gap region is the following two variation types: there is a gap, one base, Indel; there is a gap, more than one base, Indel+SNP;

[0030] The adjustment process includes: extracting an adjustment sequence, extracting an initial reference sequence, performing an adjustment comparison operation according to the adjustment sequence and the initial reference sequence, obtaining an adjustment comparison result after the adjustment comparison of the adjustment sequence, performing result embedding processing on the adjustment comparison result and the initial comparison result, and obtaining the comparison result after the result embedding processing as the final comparison result;

[0031] Optionally, the extraction of the adjustment sequence includes: boundary determination and sequence extraction;

[0032] Optionally, the adjustment and comparison operation includes: comparing the adjusted sequence with the initial reference sequence, wherein the comparison method adopts a gap-filling comparison and a gap-adding comparison;

[0033] Optionally, the result embedding process includes: deleting the old sequence between the start point and the end point determined when the adjustment sequence is cut, and embedding the adjustment comparison result between the start point and the end point determined when the adjustment sequence is cut.

[0034] The multiple sequence alignment method further comprises: integrating the final alignment result with the bridge sequence information of the pan-genome sequence data to obtain a multiple sequence alignment result of the pan-genome sequence data;

[0035] Optionally, the bridge is a single branch. If any point or edge of the single branch is deleted, the number of connected blocks in the colored graph will increase.

[0036] An analysis device for a graph-based multiple sequence alignment method, the device comprising: a memory and a processor;

[0037] The memory is used to store program instructions;

[0038] The processor is used to call program instructions, and when the program instructions are executed, it is used to execute the above-mentioned graph-based multiple sequence alignment method.

[0039] An analysis system for a graph-based multiple sequence alignment method, comprising:

[0040] A first processing unit is used to obtain a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome;

[0041] A second processing unit is used to mark and obtain the characteristics of the access state of a single node of the colored graph, and traverse the colored graph to obtain a cSupB data model after decomposing the colored graph;

[0042] The third processing unit is used to obtain the cSupB data model after the colored graph is decomposed, and extract the offset value data features of the nodes in the cSupB data model;

[0043] The fourth processing unit performs a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model.

[0044] Optionally, the analysis system also includes: an adjustment judgment unit, which is used to analyze the initial comparison result to obtain a data analysis result of the initial comparison result, compare the relationship between the data analysis result and the adjustment standard, and determine whether to adjust the initial comparison result according to the relationship, and output the final comparison result after or without adjustment; after adjusting the initial comparison result, output the adjusted comparison result as the final comparison result; without adjusting the initial comparison result, output the initial comparison result as the final comparison result; the adjustment judgment unit is located after the fourth processing unit.

[0045] Optionally, the analysis system also includes: a screening processing unit, used to screen the cSupB data model after the colored graph decomposition to identify the target cSupB data model; perform a first preprocessing on the target cSupB data model to obtain an initial comparison result; the initial comparison result is an initial cSupB comparison matrix; the screening processing unit is located after the third processing unit and before the fourth processing unit.

[0046] Optionally, the analysis system further comprises: a fifth processing unit, configured to integrate the final alignment result with the bridge sequence information of the pan-genome sequence data to obtain a multiple sequence alignment result of the pan-genome sequence data. The fifth processing unit is located after the adjustment and judgment unit.

[0047] A computer-readable storage medium stores a computer program, which implements the above-mentioned graph-based multiple sequence alignment method when executed by a processor.

[0048] This application has the following beneficial effects:

[0049] 1. The present invention proposes a research idea based on pan-genome, which depicts a large amount of complex data by constructing a colored graph with a clear data structure, and then decomposes and reorganizes the colored graph to effectively analyze the structural information of the colored graph; a set of pan-genome sequence data can be constructed into a directed acyclic graph with unique start and end nodes. Based on this graph, on the one hand, a cSupB data structure is proposed, which inherits the superbubble characteristics while also combining sample information such as node source and linkage. By traversing the graph once through the proposed post-order traversal strategy, the entire colored graph can be decomposed into large and small cSupBs. On the other hand, in order to describe the position information of the node, an offset value is introduced.

[0050] 2. The present invention proposes a multiple sequence alignment method GMSA based on a colored graph, which takes cSupB as the basic unit. On the basis of cSupB, the initial alignment result can be obtained by simply traversing the colored graph using the offset value obtained in the graph decomposition process (the preliminary alignment result is obtained by using the property of the longest path of the offset value), and then a micro-alignment operation is performed on the position where the alignment result may have errors, so as to obtain the final alignment result simply, quickly and accurately; in this process, not only the concept of the guide tree is abandoned and the length and number of the alignment sequences are reduced, but also the scoring matrix is ​​simplified and the objective function is optimized, the accuracy of the alignment is improved, and the time and space complexity of the alignment are reduced. Finally, all cSupB alignment situations and bridge sequence information are integrated to obtain the multiple sequence alignment results at the whole genome level.

[0051] The technical solution of the present application effectively solves the problem that the traditional multi-sequence alignment method cannot meet the current situation when faced with a large number of sequenced sample genomes, and the results are fast and accurate. BRIEF DESCRIPTION OF THE DRAWINGS

[0052] In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the drawings required for use in the description of the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative work.

[0053] Figure 1 A schematic flow chart of a graph-based multiple sequence alignment method provided in an embodiment of the present invention;

[0054] Figure 2 A schematic diagram of an analysis device for a graph-based multiple sequence alignment method provided by an embodiment of the present invention;

[0055] Figure 3 A schematic flow chart of an analysis system for a graph-based multiple sequence alignment method provided by an embodiment of the present invention;

[0056] Figure 4 A schematic diagram of the construction, decomposition and reorganization of a colored graph provided in an embodiment of the present invention;

[0057] Figure 5 A schematic diagram of the GMSA analysis process of a colored graph provided by an embodiment of the present invention;

[0058] Figure 6 The filling process of the initial comparison matrix provided by the embodiment of the present invention;

[0059] Figure 7 Two comparison schemes, GPA and GJA, are provided in the embodiment of the present invention; DETAILED DESCRIPTION

[0060] In order to enable those skilled in the art to better understand the solutions of the present invention, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the accompanying drawings in the embodiments of the present invention.

[0061] In some of the processes described in the specification and claims of the present invention and the above-mentioned figures, multiple operations that appear in a specific order are included, but it should be clearly understood that these operations may not be executed in the order in which they appear in this article or executed in parallel. The serial numbers of the operations, such as 101, 102, etc., are only used to distinguish different operations, and the serial numbers themselves do not represent any execution order. In addition, these processes may include more or fewer operations, and these operations may be executed in sequence or in parallel. It should be noted that the descriptions of "first", "second", etc. in this article are used to distinguish different messages, devices, modules, etc., do not represent the order of precedence, and do not limit the "first" and "second" to be different types.

[0062] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative work are within the scope of protection of the present invention.

[0063] Figure 1 A schematic flow chart of a graph-based multiple sequence alignment method provided in an embodiment of the present invention, specifically, the method comprises the following steps:

[0064] 101: Obtain a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome;

[0065] In one embodiment, after obtaining the pan-genome sequence data, preprocessing is performed, and the preprocessing adopts a method including base replacement and addition of sequence fragments;

[0066] The base replacement method is to replace the base with the base with the highest frequency in other series of the site when there is a degenerate base sequence in the pan-genome sequence data;

[0067] The added sequence fragment is to add the same sequence fragment to the head and tail of each sequence in the pan-genome sequence data; the role of this fragment is to anchor the head and tail of all sequences together, and should not affect the topological structure of other positions on the graph, that is, it will not form a ring structure with the genome sequence. Therefore, when the fragment is randomly generated, it needs to be judged whether it is applicable to this group of samples. It is not fixed and will be adjusted accordingly as the samples change.

[0068] After the preprocessing, the pan-genome colored graph is a directed, acyclic, non-degenerate base graph with a unique starting point and end point. Each sequence corresponds to a unique path on the colored graph starting from the starting point and ending at the end point.

[0069] Pan-genome is a collection of multiple genomes, including core genes and non-essential genes.

[0070] Path: from v0 to v k A path is a sequence v0, e1, v1, e2, ..., e k , v k , where e i is the connecting node v i-1 to v i The length of the path is k. If there is a path in the graph with the same starting and ending points, then the path is "closed", indicating that the graph has a cycle.

[0071] 102: Mark and obtain the feature of the access state of a single node of the colored graph, traverse the colored graph to obtain a cSupB data model after decomposing the colored graph;

[0072] In one embodiment, the colored graph is traversed to obtain data information of the cSupB data model after decomposing the colored graph; the characteristics of the access state of a single node of the colored graph include an unvisited state, a semi-visited state, an accessible state and a visited state;

[0073] The unvisited state is that no entry point has been visited, recorded as -1; the semi-visited state is that at least one entry point has been visited and at least one entry point has not been visited, recorded as 0; the accessible state is that all entry points have been visited and the node itself is in a state that can be visited, recorded as 1; the visited state is that all entry points of the node have been visited and the node itself has also been visited, recorded as 2;

[0074] The traversal of the colored graph adopts a post-order traversal method; the basic requirement of post-order traversal is that in a graph, all child nodes are visited before the parent node. The basic requirement of post-order traversal is that all parent nodes are visited before the child nodes, that is, for each node in the graph, only when all its entry points have been visited can it visit itself and its exit points.

[0075] The cSupB data model is: in a colored graph G = (V, E, C), V(G), V(E) and V(C) are the vertex set, edge set and color set of graph G respectively. For any color set C1(C1C), G1 = (V1, E1, C1) is a subgraph of graph G, satisfying that for any node u i ∈V1, For two different points s and t, it is called<s,t,C1> A colored SuperBubble, also known as cSupB; if it satisfies the four traditional superbubble criteria on G1, reachability, matching, acyclicity and minimality, the function color is used to represent the color information of the node or edge; s is called the source node, and t is the sink node; given any two cSupBs, defined as cSupB1 and cSupB2, let G1 = (V1, E1, C1) and G2 = (V2, E2, C2) be the subgraphs induced by the nodes contained in cSupB1 and cSupB2 respectively. Then let V0 = V1∩V2, C0 = C1∩C2, when or and When , the outputs cSupB1 and cSupB2 are separated; when and When V0≠V1(V2), the output cSupB1 and cSupB2 intersect; when V0=V1 and C0=C1, or V0=V2 and C0=C2, the output cSupB1 and cSupB2 are included. If two cSupBs are in an inclusion relationship, for example: V0=V2 and C0=C2, cSupB1 is called the parent bubble of cSupB2, and cSupB2 is the child bubble of cSupB1.

[0076] The data model of SupB is: in the assembly graph G = (V, E), V(G) and V(E) are the point set and edge set of graph G respectively. For any two different points s and t,<s,t,> It is called a SuperBubble; if the above SuperBubble satisfies the following four criteria: reachability, matching, acyclicity and minimality; reachability is denoted as reachability: there is a path from point to point; matching is denoted as matching: the set of points that can be reached from point without passing through point is the same as the set of points that can reach point without passing through point; acyclicity is denoted as acyclicity: the subgraph induced by the point set that satisfies matching is acyclic; minimality is denoted as minimality: in the point set U, except point t, no other point can form a pair with point s that can satisfy the above three criteria; s is called the source node, and t is the sink node; here, only superbubbles containing at least two supernodes are considered.

[0077] like Figure 4 As shown, Figure 4 Ahaplotype is a single type specimen, which is three initial sequences; based on the SupB data model, Figure 4 There are only two superbubbles in B:<TAT,ACC> and<TCA,GTA> ; Figure 4 B is a colored graph, k=3, each circle represents a node, and the black arrow represents an edge. The characters above the nodes represent the base, access order, and color of the nodes in turn, and the numbers at the bottom represent the theoretical offset value of each node; Figure 4 C is the result of accessing the graph in the reverse direction; the numbers under the nodes represent the final and initial offset values ​​respectively; Figure 4 D is the cSupB structure tree. Five cSupBs are found here, followed by bub1.<TAT,ACC,111> ;bub2.<GGG,GTA,110> ;bub3.<CAC,GTA,011> ;bub4.<TCA,GGG,110> ; and bub5.<TCA,GTA,111> ; Figure 4 E is the final decomposition and representation of the colored graph.

[0078] The data information of the cSupB data model includes but is not limited to the following information: source point, sink point, color of the cSupB data model and order of the cSupB data model.

[0079] The matching principle of the cSupB data model: the colors of adjacent outgoing (incoming) edges of the source (sink) node do not intersect; the colors of adjacent incoming (outgoing) edges of the sink (source) node and the intersection of the colors of the outgoing (incoming) edges of each adjacent source (sink) node are not empty; the cSupB color is the intersection of the colors of the source node and the sink node.

[0080] Matching process: In the post-order traversal of the class, if a source node s is encountered, it is placed in a queue of source nodes to be visited Q; if a sink node t is encountered, it starts to find its matching source node s from Q in reverse according to the matching principle. In particular, if it is Then remove s from Q and continue; if it is Then stop matching.

[0081] Incoming node / incoming edge: Given any two nodes u and v, u is called an incoming node of v if and only if there is at least one path from u to v, and all the edges on the path are called incoming edges of v. The number of adjacent incoming edges of v is called the indegree of v.

[0082] Outgoing node / outgoing edge: At the same time, v is also called the outgoing node of u, all the edges on the path are called the outgoing edges of u, and the number of adjacent outgoing edges of u is called the outdegree of u.

[0083] Degree: The sum of the in-degree and out-degree of a node is called the degree of the node.

[0084] supernode: If at least one of the out-degree or in-degree of a node is greater than 1, then the node is called a supernode.

[0085] branch: In this article, a path v0, e1, v1, e2, ..., e k , v k It is called a branch if the out-degree and in-degree of all nodes on the path are equal to 1 and the previous node of v0 and v k The next node after is a supernode.

[0086] Bridge: In this article, a branch is called a bridge. If any vertex or edge on the branch is deleted, the number of connected blocks in the graph will increase.

[0087] Bubble: Due to sequence differences, the pathways on the graph diverge first and then converge, forming a bubble-like structure called a bubble. For details, see the SupB data model.

[0088] 103: Obtain the cSupB data model after the colored graph is decomposed, and extract the offset value data features of the nodes in the cSupB data model.

[0089] In one embodiment, the offset value data refers to the longest path of the offset value; for any directed acyclic graph G = (V, E), if for any two nodes u, v∈V and there is at least one path between the two points, let pos(u) = a, pos(v) = b, then |ab| is the maximum distance between u and v.

[0090] 104: Perform a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model.

[0091] In one embodiment, the first preprocessing includes: defining the cSupB data model as an s×l cSupB alignment matrix based on the offset value data feature, traversing the cSupB alignment matrix, adding the base information of each node to the corresponding position of the cSupB alignment matrix, and obtaining an initial alignment result of the cSupB data model;

[0092] Optionally, in the cSupB comparison matrix, the cSupB comparison matrix is ​​denoted as M, s is the number of samples of the target cSupB data model, l=pos(source)-pos(sink)+1, and l is the maximum length of the target cSupB data model;

[0093] Optionally, the traversal process includes but is not limited to the following traversal methods: a post-order traversal method; starting from the source point, traversing along the direction of the edge of the graph G, and adding the base information of each node (the last base of the k-mer) to the corresponding position of M.

[0094] Screening the cSupB data model to identify a target cSupB data model;

[0095] Optionally, the screening process includes: extracting a root cSupB data model of level=1 in the cSupB data model; the root cSupB data model is the target cSupB data model;

[0096] Optionally, the target cSupB data model is preprocessed for the first time to obtain an initial comparison result.

[0097] like Figure 6 As shown, A is the initial colored graph and the selected cSupB, and B is the unfilled, filled, and filled initial alignment matrix.

[0098] 105: Analyze the initial comparison result to obtain a data analysis result of the initial comparison result, compare the relationship between the data analysis result and the adjustment standard, and determine whether to adjust the initial comparison result according to the relationship; perform or not perform adjustment processing, and output the final comparison result.

[0099] Optionally, the process of analyzing the initial alignment result includes: identifying the variation type of each aligned site in the initial alignment result, analyzing the initial alignment result column by column to obtain a hash table pos2vartype, traversing the hash table pos2vartype in sequence according to the position size, and obtaining a data analysis result;

[0100] Optionally, the variation types include the following five types: no gap, 1 base, no variation; no gap, 2 bases, SNP; no gap, 3 or 4 bases, SNP; with gap, 1 base, Indel; with gap, more than 1 base, Indel+SNP;

[0101] Optionally, the key of the hash table pos2vartype is the matrix position, and the value is the variation type;

[0102] Optionally, the method for traversing the hash table pos2vartype is: traverse in sequence according to the size of the matrix position, encounter a gap area during the traversal, the output length is the gap length, and the unaligned score of the m bp sequence (Exploring Range; ER, default: 10) on both sides of the gap area is calculated, and the unaligned score (UnAligned Score; UAS) is the data analysis result;

[0103]

[0104] if

[0105] It means that there is a high possibility of comparison deviation in the gap area, that is, the data analysis results meet the adjustment standard; g0 and m / 2 are called adjustment standards; the standard can be adjusted by changing the value. The default parameters in this article are g0=2 and m=10.

[0106] If the gap region does not need fine-tuning, there is no gap region or there is a gap region but it does not meet the fine-tuning standard, that is, no adjustment is required, and the current alignment matrix M is the final alignment result.

[0107] Optionally, the variation type of the gap region is the following two variation types: there is a gap, one base, Indel; there is a gap, more than one base, Indel+SNP.

[0108] The adjustment process includes: extracting an adjustment sequence, extracting an initial reference sequence, performing an adjustment comparison operation according to the adjustment sequence and the initial reference sequence, obtaining an adjustment comparison result after the adjustment comparison of the adjustment sequence, performing result embedding processing on the adjustment comparison result and the initial comparison result, and obtaining the comparison result after the result embedding processing as the final comparison result;

[0109] Optionally, the extraction of the adjustment sequence includes: boundary determination and sequence extraction; boundary determination: starting from the gap region selected in the previous step, start to extend to both sides until there are continuous m0 (Anchor Length; AL, default: 10) positions without variation, and their farthest ends are respectively used as the start and end positions of the fine-tuning region, that is, the boundary of the fine-tuning region is determined. If the boundary of the previous fine-tuning is encountered during the extension process, stop the extension and use the previous position as the boundary point of the new fine-tuning region, that is, ensure that different fine-tuning regions do not intersect with each other.

[0110] Sequence extraction: extract all sequences in the region from the initial alignment results, perform deduplication operations to obtain a set of unique sequences for subsequent alignment. Let's assume that the unique sequences obtained are s0 (s0 < s), and save these sequences. In addition to saving the sequences, it is also necessary to save information such as the samples contained and the number of gaps. The way to save is to use a hash table with the extracted sequence as the key and the sample information or the number of gaps as the value.

[0111] Determination of the initial reference sequence: Starting from the source node, a path with continuous offset values ​​is randomly selected. This path must exist and its length is pos(sink)-pos(source)+1. It is recorded as the initial reference sequence. This sequence must be the longest sequence in the region at this time. The length of other sequences must be less than or equal to the reference sequence. Of course, this reference sequence may not exist in the sample.

[0112] The extracted sequences are sorted in a priority manner from small to large gap numbers and from many to few sample numbers, and gaps in the sorted sequences are deleted; the sorting process is actually an evaluation process of the importance of the sequences, and the sorting result is actually also an evaluation result of the importance of the sequences. In addition, in the subsequent micro-alignment process, the score matrix can be dynamically adjusted according to the importance of the sequences.

[0113] Optionally, the adjustment and alignment operation includes: aligning the adjusted sequence with the initial reference sequence, using gap-filling alignment and gap-joining alignment as the alignment method; aligning the arranged sequence s' with the reference sequence s0 in turn. If the length of the sequence s' is less than s0, adopt the GPA alignment scheme to add a gap to s' so that the two sequences are of equal length; if the lengths are the same, compare the similarities between sequences s' and s0. If the similarity is high, no alignment will be performed. Otherwise, adopt the GJA alignment scheme to consider whether adding a gap will produce a better alignment result. If a new gap is added to the alignment result, it means that the longest reference sequence s0 has changed. The reference sequence s0 is updated and the fine-tuning operation is performed again from the beginning. There are two micro-alignment schemes involved here, gap-polishing alignment (GPA) and gap-joining alignment (GJA), which are applied to different alignment situations respectively. Figure 7 As shown in the figure, the application scenarios of the two comparison schemes GPA and GJA are shown respectively.

[0114] Optionally, the result embedding process includes: deleting the old sequence between the start point and the end point determined when the adjustment sequence is cut, and embedding the adjustment alignment result between the start point and the end point determined when the adjustment sequence is cut. The distribution of the five types of variation in the final alignment result is counted again, and compared with the statistical results before the adjustment process and discussed.

[0115] 106: Integrate the final alignment result with the bridge sequence information of the pan-genome sequence data to obtain a multiple sequence alignment result of the pan-genome sequence data;

[0116] Optionally, the bridge is a single branch. If any point or edge of the single branch is deleted, the number of connected blocks in the colored graph will increase.

[0117] like Figure 5 The figure shows the schematic diagram of the GMSA analysis process, which includes the steps shown in the figure. The basic alignment idea of ​​GMSA is: nesting and decomposing the directed acyclic colored graph constructed by a set of sequences to obtain the required cSupB and the initial alignment results, and then based on the initial alignment matrix, checking column by column to determine whether to fine-tune the alignment, if necessary, determine the interval and extract the required alignment sequence, and perform micro-alignment with the reference sequence (the sequence corresponding to the longest path) in turn. Finally, the obtained alignment result is re-embedded into the original alignment matrix to obtain the final alignment result and evaluate it.

[0118] Figure 2 A schematic diagram of an analysis device for a graph-based multiple sequence alignment method provided by an embodiment of the present invention, the device comprising: a memory and a processor;

[0119] The memory is used to store program instructions;

[0120] The processor is used to call program instructions, and when the program instructions are executed, it is used to execute the above-mentioned graph-based multiple sequence alignment method.

[0121] Figure 3 A schematic flow chart of an analysis system of a graph-based multiple sequence alignment method provided for an embodiment of the present invention comprises: 301: a first processing unit, configured to obtain a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome;

[0122] 302: a second processing unit, configured to mark and obtain a feature of an access state of a single node of the colored graph, and traverse the colored graph to obtain a cSupB data model after decomposing the colored graph;

[0123] 303: a third processing unit, configured to obtain a cSupB data model after the colored graph is decomposed, and extract data features of offset values ​​of nodes in the cSupB data model;

[0124] 304: a fourth processing unit, performing a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model;

[0125] Optionally, the analysis system also includes: an adjustment judgment unit, which is used to analyze the initial comparison result to obtain a data analysis result of the initial comparison result, compare the relationship between the data analysis result and the adjustment standard, and determine whether to adjust the initial comparison result according to the relationship, and output the final comparison result after or without adjustment; after adjusting the initial comparison result, output the adjusted comparison result as the final comparison result; without adjusting the initial comparison result, output the initial comparison result as the final comparison result; the adjustment judgment unit is located after the fourth processing unit.

[0126] Optionally, the analysis system also includes: a screening processing unit, used to screen the cSupB data model after the colored graph decomposition to identify the target cSupB data model; perform a first preprocessing on the target cSupB data model to obtain an initial comparison result; the initial comparison result is an initial cSupB comparison matrix; the screening processing unit is located after the third processing unit and before the fourth processing unit.

[0127] Optionally, the analysis system further comprises: a fifth processing unit, configured to integrate the final alignment result with the bridge sequence information of the pan-genome sequence data to obtain a multiple sequence alignment result of the pan-genome sequence data. The fifth processing unit is located after the adjustment and judgment unit.

[0128] A computer-readable storage medium stores a computer program, which implements the above-mentioned graph-based multiple sequence alignment method when executed by a processor.

[0129] The validation results of this validation example show that assigning intrinsic weights to indications can moderately improve the performance of the method relative to the default setting.

[0130] Those skilled in the art can clearly understand that, for the convenience and brevity of description, the specific working processes of the systems, devices and units described above can refer to the corresponding processes in the aforementioned method embodiments and will not be repeated here.

[0131] In the several embodiments provided in the present application, 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 the units is only a logical function division. There may be other division methods in actual implementation, such as 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.

[0132] 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.

[0133] In addition, each functional unit in each embodiment of the present invention 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.

[0134] A person of ordinary skill in the art can understand that all or part of the steps in the various methods of the above embodiments can be completed by instructing related hardware through a program, and the program can be stored in a computer-readable storage medium, and the storage medium may include: a read-only memory (ROM), a random access memory (RAM), a disk or an optical disk, etc.

[0135] A person skilled in the art can understand that all or part of the steps in the above-mentioned embodiment method can be completed by instructing related hardware through a program, and the program can be stored in a computer-readable storage medium. The above-mentioned storage medium can be a read-only memory, a disk or an optical disk, etc.

[0136] The above is a detailed introduction to a computer device provided by the present invention. For a person skilled in the art, according to the concept of the embodiments of the present invention, there may be changes in the specific implementation method and application scope. In summary, the content of this specification should not be understood as limiting the present invention.

Claims

1. A graph-based multiple sequence alignment method comprising: Acquire a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome; Mark and obtain the characteristics of the access status of a single node of the colored graph, traverse the colored graph to obtain a cSupB data model after decomposing the colored graph; The cSupB data model is: G=(V, E, C) middle, V(G) , V(E) and V(C) They are respectively G The point set, edge set and color set of C 1 ( C 1 C ), G 1 =(V 1 , E 1 , C 1 ) It is a picture G A subgraph of , satisfying that for any node , ; for two different points s and t , called < s , t, C 1 > a colored SuperBubble; s is called the source node, and t is the sink node; Obtaining a cSupB data model after the colored graph is decomposed, and extracting offset value data features of nodes in the cSupB data model; the offset value data refers to the longest path of the offset value; Performing a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model; The first preprocessing includes: defining the cSupB data model as an s×l cSupB comparison matrix based on the offset value data feature, traversing the cSupB comparison matrix, adding the base information of each node to the corresponding position of the cSupB comparison matrix, and obtaining an initial comparison result of the cSupB data model.

2. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The characteristics of the access state of a single node of the colored graph include an unvisited state, a semi-visited state, an accessible state and a visited state.

3. The graph-based multiple sequence alignment method according to claim 1, characterized in that: In the cSupB comparison matrix, s is the number of samples of the target cSupB data model, l=pos(source)-pos(sink)+1, and l is the maximum length of the target cSupB data model.

4. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The traversal process includes the following traversal methods: a class post-order traversal method.

5. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The base information of each node is the last base of the k-mer.

6. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The cSupB data model is screened to identify a target cSupB data model.

7. The graph-based multiple sequence alignment method according to claim 6, characterized in that: The screening process includes: extracting a root cSupB data model with level=1 in the cSupB data model; the root cSupB data model is the target cSupB data model.

8. The graph-based multiple sequence alignment method according to claim 6, characterized in that: The target cSupB data model is preprocessed for the first time to obtain an initial comparison result.

9. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The multiple sequence alignment method further comprises: analyzing the initial alignment result to obtain a data analysis result of the initial alignment result, comparing the relationship between the data analysis result and an adjustment standard, determining whether to perform adjustment processing on the initial alignment result according to the relationship; performing or not performing adjustment processing, and outputting a final alignment result.

10. The graph-based multiple sequence alignment method according to claim 1, characterized in that: The process of analyzing the initial alignment result includes: identifying the variation type of each aligned site in the initial alignment result, analyzing the initial alignment result column by column to obtain a hash table pos2vartype, traversing the hash table pos2vartype, and obtaining a data analysis result.

11. The graph-based multiple sequence alignment method according to claim 10, characterized in that: The variation types include: no gap, 1 base, no variation; no gap, 2 bases, SNP; no gap, 3 or 4 bases, SNP; gap, 1 base, Indel; gap, more than 1 base, Indel+SNP.

12. The graph-based multiple sequence alignment method according to claim 10, characterized in that: The key of the hash table pos2vartype is the matrix position, and the value is the variation type.

13. The graph-based multiple sequence alignment method according to claim 10, characterized in that: The method for traversing the hash table pos2vartype is: traversing in sequence according to the size of the matrix position, encountering a gap area during the traversal, the output length is the gap length, calculating the misalignment score of the m bp sequences on both sides of the gap area, and the misalignment score is the data analysis result.

14. The graph-based multiple sequence alignment method according to claim 13, characterized in that: The variation types of the gap region are the following two variation types: there is a gap, one base, Indel; there is a gap, more than one base, Indel+SNP.

15. The graph-based multiple sequence alignment method according to claim 9, characterized in that: The adjustment process includes: extracting an adjustment sequence, extracting an initial reference sequence, performing an adjustment comparison operation according to the adjustment sequence and the initial reference sequence, obtaining an adjustment comparison result after adjusting and comparing the adjustment sequence, performing result embedding processing on the adjustment comparison result and the initial comparison result, and obtaining the comparison result after the result embedding processing as the final comparison result.

16. The graph-based multiple sequence alignment method according to claim 15, characterized in that: The extracting and adjusting sequence includes: boundary determination and sequence extraction.

17. The graph-based multiple sequence alignment method according to claim 15, characterized in that: The adjustment and comparison operation includes: comparing the adjusted sequence with the initial reference sequence, and the comparison method adopts the gap filling comparison and gap adding comparison.

18. The graph-based multiple sequence alignment method according to claim 15, characterized in that: The result embedding process includes: deleting the old sequence between the start point and the end point determined when the adjustment sequence is cut, and embedding the adjustment comparison result between the start point and the end point determined when the adjustment sequence is cut.

19. The graph-based multiple sequence alignment method according to claim 15, characterized in that: The multiple sequence alignment method further includes: integrating the final alignment result with the bridge sequence information of the pan-genome sequence data to obtain a multiple sequence alignment result of the pan-genome sequence data.

20. The graph-based multiple sequence alignment method according to claim 19, characterized in that: The bridge is a single branch. If any point or edge of the single branch is deleted, the number of connected blocks in the colored graph will increase.

21. An analysis device for a graph-based multiple sequence alignment method, the device comprising: Memory and processor; The memory is used to store program instructions; The processor is used to call program instructions, and when the program instructions are executed, it is used to execute the graph-based multiple sequence alignment method described in any one of claims 1-20.

22. An analysis system for a graph-based multiple sequence alignment method, comprising: A first processing unit is used to obtain a set of pan-genome sequence data, map the pan-genome sequence data, and obtain a colored map of the pan-genome; A second processing unit is used to mark and obtain the characteristics of the access state of a single node of the colored graph, and traverse the colored graph to obtain a cSupB data model after decomposing the colored graph; The cSupB data model is: G=(V, E, C) middle, V(G) , V(E) and V(C) They are respectively G The point set, edge set and color set of C 1 ( C 1 C ), G 1 =(V 1 , E 1 , C 1 ) It is a picture G A subgraph of , satisfying that for any node , ; for two different points s and t , called < s , t, C 1 > a colored SuperBubble; s is called the source node, and t is the sink node; The third processing unit is used to obtain the cSupB data model after the colored graph is decomposed, and extract the offset value data features of the nodes in the cSupB data model; the offset value data refers to the longest path of the offset value; A fourth processing unit, performing a first preprocessing on the cSupB data model based on the offset value data feature to obtain an initial comparison result of the cSupB data model; The first preprocessing includes: defining the cSupB data model as an s×l cSupB comparison matrix based on the offset value data feature, traversing the cSupB comparison matrix, adding the base information of each node to the corresponding position of the cSupB comparison matrix, and obtaining an initial comparison result of the cSupB data model.

23. A computer-readable storage medium having a computer program stored thereon, wherein when the computer program is executed by a processor, the graph-based multiple sequence alignment method according to any one of claims 1 to 20 is implemented.

Images