Error correction code table generation method and encoding and decoding method based on generated code table

Abstract

The invention provides an error correction code table generation method and a coding and decoding method based on a generated code table. The error correction code table generation method comprises the following steps: determining a limited character table of a sequence code serving as an error correction target and a code table size of a corresponding error correction code, and generating an error correction feature sequence based on the limited character table and an error correction error type, the error correction feature sequence is used for defining the number of error correction errors of the code table based on the types of the error correction errors; determining a feature sequence, a feature mask and a translation amount based on the error correction feature sequence; iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount; and in response to the condition that the loss function of the dynamic programming matrix is smaller than or equal to a preset threshold value, determining two code words in a code table corresponding to each dynamic programming matrix subjected to iteration completion, so as to obtain an error correction code table. In this way, high-performance error correction code table generation can be realized.

Description

Technical Field

[0001] The present application relates to the technical field of data error correction, and more specifically, to a method for generating an error correction code table and a coding and decoding method based on the generated code table. Background Art

[0002] Sequence code storage, such as DNA storage, is an emerging technology that increasingly pursues high density, long-term stability and large capacity. In the process of sequence code storage, complex errors occur due to the storage environment, making the design of error correction codes challenging.

[0003] Here, the basic idea of ​​error correction code can be understood as defining a set of (one or more) functions in the sequence space, through which the received sequence is evaluated. These functions are usually used to detect and correct errors in the sequence. If the function value of the received sequence meets expectations, that is, matches the pre-defined rules, then it can be regarded as a correct sequence; otherwise, judgment or correction is made based on different function value results to restore the original information as much as possible.

[0004] Therefore, it is desirable to provide an improved error correction code table generation technique to generate an error correction code table with high performance. Summary of the invention

[0005] An embodiment of the present application provides a method for generating an error correction code table, which iterates a dynamic programming matrix based on an error correction feature sequence that defines the number of error correction errors of a predetermined error correction error type so that any two codewords of the code table corresponding to the dynamic programming matrix can be distinguished based on the number of error correction errors of the predetermined error correction error type, thereby achieving high-performance error correction code table generation.

[0006] The embodiment of the present application further provides a coding and decoding method based on generating a code table, which realizes high-performance sequence coding and decoding through the error correction capability of the code table based on predetermined error correction error types.

[0007] According to one aspect of the present application, a method for generating an error correction code table is provided, comprising: determining a finite character table of a sequence code as an error correction target and a code table size of a corresponding error correction code, and generating an error correction feature sequence based on the finite character table and the error correction error type, the error correction feature sequence being used to define the number of error correction errors of the code table based on the error correction error type; determining a feature sequence, a feature mask and a shift amount based on the error correction feature sequence; iterating a dynamic programming matrix based on the feature sequence, the feature mask and the shift amount; and, in response to a loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, determining two codewords in the code table corresponding to each dynamic programming matrix completed by the iteration to obtain an error correction code table.

[0008] According to another aspect of the present application, a coding method based on a generated code table is provided, comprising: determining a floor value of a logarithmic function with base 2 of a code table size of the generated code table; preprocessing an input binary sequence to perform zero padding at the end and recording the number of zero paddings; and group encoding the preprocessed input binary sequence using the number of the floored values ​​as a group.

[0009] According to another aspect of the present application, a decoding method based on a generated code table is provided, comprising: dividing a received sequence to be decoded into a predetermined number of groups, the predetermined number being the number of coding groups, and the length of each group being the length of each codeword in the generated code table; comparing each group with the code table starting from the first group to determine the codeword corresponding to the group in the code table; for the kth group in the sequence to be decoded and the corresponding error correction sequence index being i, adding a prefix of the codeword minus a window size to the codeword in the generated code table, and matching the codeword with a sequence of the length of the previous correction sequence minus the start of the window size and the length of the codeword multiplied by k; and, when all groups of the sequence to be decoded are matched, concatenating all correct sequences to obtain a decoding sequence.

[0010] According to another aspect of the present application, a method for generating an error correction code table for a gene sequence is provided, comprising: determining the size of the error correction code table for the gene sequence, and generating an error correction feature sequence based on the gene sequence {A, T, C, G} as a finite character table and the error correction error types and the number of error correction errors corresponding to the gene sequence error correction target; determining a feature sequence, a feature mask and a translation amount based on the error correction feature sequence; iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount; and, in response to a loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, determining two codewords in the code table corresponding to each dynamic programming matrix completed by the iteration, so as to obtain an error correction code table for the gene sequence.

[0011] According to another aspect of the present application, a coding method based on a gene sequence code table is provided, comprising: determining a floor value of a logarithmic function with base 2 of the code table size of the gene sequence code table; preprocessing an input binary sequence to pad the end with zeros and recording the number of padding zeros; and grouping and encoding the preprocessed input binary sequence according to the number of the floor values ​​as a group.

[0012] According to another aspect of the present application, a decoding method based on a gene sequence code table is provided, comprising: dividing a received gene sequence to be decoded into a predetermined number of groups, the predetermined number being the number of coding groups, and the length of each group being the length of each codeword in the gene sequence code table; comparing each group with the gene sequence code table starting from the first group to determine the codeword corresponding to the group in the gene sequence code table; for the kth group in the gene sequence to be decoded and the corresponding error correction sequence index being i, adding a prefix of the codeword minus a window size to the codeword in the gene sequence code table, and matching the codeword with a sequence of the length of the previous correction sequence minus the start of the window size and the length of the codeword multiplied by k; and, when all groups of the gene sequence to be decoded are matched, splicing all correct sequences to obtain a decoded gene sequence.

[0013] The error correction code table generation method provided in the embodiment of the present application can iterate the dynamic programming matrix based on the error correction feature sequence that defines the number of error correction errors of a predetermined error correction error type so that any two codewords of the code table corresponding to the dynamic programming matrix can be distinguished based on the number of error correction errors of the predetermined error correction error type, thereby realizing high-performance error correction code table generation.

[0014] Furthermore, the encoding and decoding method based on the gene sequence code table provided in the embodiment of the present application can achieve high-performance gene sequence encoding and decoding through the error correction capability of the gene sequence code table based on predetermined error correction types. BRIEF DESCRIPTION OF THE DRAWINGS

[0015] By reading the detailed description of the preferred specific embodiments below, various other advantages and benefits of the present application will become clear to those of ordinary skill in the art. The drawings in the specification are only used for the purpose of illustrating the preferred embodiments and are not considered to be limitations of the present application. Obviously, the drawings described below are only some embodiments of the present application. For those of ordinary skill in the art, other drawings can also be obtained based on these drawings without creative work. Moreover, the same reference numerals are used to represent the same components throughout the drawings.

[0016] Figure 1 The figure illustrates a schematic flow chart of a method for generating an error correction code table according to an embodiment of the present application.

[0017] Figure 2 The figure illustrates a schematic flow chart of a coding method based on a generated code table according to an embodiment of the present application.

[0018] Figure 3 The figure illustrates a schematic flow chart of a decoding method based on a generated code table according to an embodiment of the present application.

[0019] Figure 4The figure illustrates a schematic flow chart of a method for generating an error correction code table for a gene sequence according to an embodiment of the present application.

[0020] Figure 5 The figure illustrates a schematic flow chart of an encoding method based on a gene sequence code table according to an embodiment of the present application.

[0021] Figure 6 The figure illustrates a schematic flow chart of a decoding method based on a gene sequence code table according to an embodiment of the present application.

[0022] Figure 7 A schematic diagram of an application example of a method for generating an error correction code table for a gene sequence and a coding and decoding method based on a gene sequence code table according to an embodiment of the present application is illustrated.

[0023] Figure 8 The diagram illustrates an application example of the method for generating an error correction code table for a gene sequence according to an embodiment of the present application.

[0024] Fig. 9 The diagram illustrates an application example of the encoding method based on the gene sequence code table according to an embodiment of the present application.

[0025] Fig.10 The diagram illustrates an application example of a decoding method based on a gene sequence code table according to an embodiment of the present application.

[0026] Fig.11 The diagram illustrates an application example of block-by-block comparison in a decoding method based on a gene sequence code table according to an embodiment of the present application. DETAILED DESCRIPTION

[0027] Below, the exemplary embodiments according to the present application will be described in detail with reference to the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments of the present application, and it should be understood that the present application is not limited to the exemplary embodiments described here.

[0028] Schematic code table generation method

[0029] Figure 1 FIG. 1 is a schematic flow chart of a method for generating an error correction code table according to an embodiment of the present application. Figure 1 As shown, the method for generating an error correction code table according to an embodiment of the present application includes the following steps.

[0030] Step S110, determine the finite alphabet of the sequence code as the error correction target and the code table size of the corresponding error correction code, and generate an error correction feature sequence based on the finite alphabet and the error correction error type, wherein the error correction feature sequence is used to define the number of error correction errors of the code table based on the error correction error type.

[0031] Here, the error correction code table generation method according to the embodiment of the present application can be applied to any form of sequence code, such as a binary sequence {0,1}, a p-binary sequence {0,1,…,p-1}, a genetic sequence {A,T,C,G}. This sequence code can be represented as a finite character table, for example, denoted as Σ, and where |Σ|<+∞, which is a set composed of a finite number of different character elements, as long as different character elements can be distinguished from each other, so as to indicate the characters of the code used.

[0032] In addition, since the code table can be regarded as a two-dimensional character array, the number of code words in the code table is n, and the length of each code word is l, then the size of the code table is n×l. Here, in the field of coding technology, there are two basic principles. First, the code table has a certain error correction capability, which is equivalent to a sufficiently large distance between any two code words in the code table; second, as long as different suffixes are continuously added between any two code words, the distance between them can be arbitrarily large.

[0033] Therefore, compared with the classic coding idea, which is to design a certain specific error correction rule to equivalently obtain a code table that meets the corresponding requirements, the applicant of this application considers that since the focus is on the code table, if the design of the rule is very difficult, the code table can be directly constructed to implement specific error correction rules based on the above principles.

[0034] Furthermore, considering the considerations based on different error correction requirements and the principle of implementing error correction by adding different suffixes between any two code words, the applicant of the present application generates a code table based on a "bottom-up" growth method, that is, in the error correction code table generation method according to the embodiment of the present application, after determining the predetermined code table size, it is continuously grown until the error correction requirements are met. Here, the error correction requirements may include the types of errors that need to be corrected, and the number of errors that need to be corrected for each error.

[0035] In the embodiment of the present application, a complete character table can be defined based on the limited character table, which is expressed as Among them, * represents a wildcard. When used for character matching, it is equal to any element in Σ by default; = means that the characters of the two sequences corresponding to the position of the symbol should be equal.

[0036] Therefore, represents non-intersection, and the feature string space can be defined based on the complete character table, which is expressed as Indicates n All permutations and combinations of elements in Represents an empty string. That is, the feature string space represents the full set of sequences composed of all finite complete character tables, and a monoid can be defined on S by character concatenation, namely:

[0037]

[0038] |p|=l represents the length of p.

[0039] In this way, based on two different elements in the monoid defined on the feature string space, a feature sequence pair representing the type of error correction can be defined, namely: (p1, p2)∈S×S=S 2 . And, for different p i ∈S corresponds to different types of error correction. For example, based on common errors in sequence codes, they can be roughly divided into the following four types:

[0040] (p1,p2),p i ≠-, i = 1, 2, represents replacement;

[0041] (-,p2),p2≠-, represents insertion;

[0042] (p1,-), p1≠-, represents deletion; and

[0043] (-,-), indicating that the mismatch is not within the error correction range considered.

[0044] As for other types of errors, such as swapping, they can also be considered as a substitution, for example: the DNA sequence AGCT becomes TGCA, which meets the judgment of the sequence pair (A==T, T==A). For another example, for the burst situation of t insertion-s deletion, substitution can usually be used to achieve it, such as: (p1, p2), to complete the match.

[0045] Then, based on the characteristic sequence pair representing the error type of error correction, an error correction characteristic unit can be defined, which is expressed as:

[0046] The first coordinate in the error correction feature unit is a sequence of t≥1 feature pairs (the actual result will not change due to the change of the sequence), which can be called a unit feature; and the second coordinate represents the maximum number of occurrences of the set of feature pairs, which can be called an allowable upper limit; Represents an integer greater than or equal to zero. If Where S 0 Represents an empty set tuple, then the error correction feature unit can be defined The "AND" operation on , that is:

[0047]

[0048] Where S * The addition of is the union.

[0049] And, for the edit distance, it can be represented by the following error correction feature unit:

[0050] (((*,*),(-,*),(*,-)),n)

[0051] Indicates that at most n edit distances are considered.

[0052] Then, based on the error correction feature unit, the error correction feature sequence is defined as: represents a sequence consisting of m (finite) error correction feature units, and the i-th unit has t i In addition, the "OR" operation can be defined on the error correction feature sequence, that is, the error correction feature units are directly combined together. For the "AND" operation, it is completed by applying the distributive law to each unit, and finally presented in the form of "conjunctive normal form".

[0053] In the embodiment of the present application, the purpose of the error correction feature sequence is to set the error correction capability of the code table, and to determine the bit array size, feature mask and translation required in the subsequent state transfer through each maximum feature number. In addition, the error correction feature sequence is also used to determine the size of the dynamic programming matrix, and to determine whether the current position meets the conditions for state transfer when a mismatch occurs through a feature sequence pair group. This will be described in detail below.

[0054] In addition, when constructing the code table, it is necessary to ensure that no two sequences in the code table can be transformed from one to the other. Therefore, when setting the error correction feature unit, if an asymmetric situation occurs, it needs to be symmetric. Specifically, the symmetric processing includes the following steps: First, for the error correction feature sequence, start traversing each error correction feature unit c i ; Then, for each feature unit c i , multiply the allowable upper limit by 2, that is, c i [1]←2*c i [1]; Finally, for the feature unit c i Every characteristic pair (p1,p2)∈c i [0], if (p2,p1)∈c i [0], then skip, otherwise, put (p2,p1) into c i [0] In.

[0055] It is worth noting that if the received sequence needs to be compared, then there is no need to consider how the received sequence is converted into a sequence in the code table, and it is only necessary to consider whether a sequence in the code table can be converted into the received sequence.

[0056] Therefore, in the error correction code table generation method according to the embodiment of the present application, generating an error correction feature sequence based on the finite character table and the error correction error type includes: generating a complete character table based on the finite character table; defining a feature string space through the non-intersecting union of multiple complete character tables, and the feature string space contains an empty character string representation; defining a feature sequence pair representing the error correction error type through two different elements in a monoid defined in the feature string space by character splicing; defining an error correction feature unit based on multiple feature sequence pairs and an allowable upper limit that is a positive integer greater than or equal to zero; and, defining an error correction feature sequence based on a sequence composed of multiple error correction feature units.

[0057] Furthermore, in the above-mentioned error correction code table generation method, the characteristic sequence pairs indicating the types of error correction errors include:

[0058] represents the replacement (p1, p2), p i ≠-,i=1,2;

[0059] represents the insertion of (-,p2), p2≠-; and

[0060] Represents deleted (p1,-), p1≠-;

[0061] Where p1, p2 are two different elements in the monoid, and - is the empty string.

[0062] In addition, in the above-mentioned error correction code table generation method, after generating the error correction feature sequence, it further includes symmetrizing the error correction feature sequence, including: starting to traverse each error correction feature unit for the error correction feature sequence; multiplying the allowable upper limit of each error correction feature unit by two; for each feature sequence pair in the error correction feature unit, determining whether a symmetric feature sequence pair that is arranged opposite to two elements of the feature sequence is in the error correction feature unit, and in response to the symmetric feature sequence pair not being in the error correction feature unit, placing the symmetric feature sequence pair in the error correction feature unit.

[0063] Step S120, determining a feature sequence, a feature mask and a translation amount based on the error correction feature sequence. The feature sequence is a sequence obtained by each error correction feature unit in the error correction feature sequence according to the first projection, that is, for any given error correction feature sequence as described above The sequence (S * 1,…,S * m ), that is, the sequence of unit features, is used as the feature sequence for comparison in subsequent dynamic programming matrix iterations.

[0064] Furthermore, when the edit number information is stored in a bit array, for any given error correction feature sequence as described above Obtain the array (n1,…,n m ), that is, according to the order of the error correction feature sequence, obtain the array about the allowable upper limit. And, define the length of the bit array as Thus, if the determination reaches the (i, j)th step (i.e., examining all possible errors that occur from the first i position of the first string to the first j position of the second string), one of the error numbers can be represented by the tuple (x1,…,x m ),0≤x i ≤n i To represent, it is called the editing information of the (i, j)th position.

[0065] Therefore, we want to use this tuple to correspond to a certain index of the bit array, where the value of the index is 1, and for indexes that are not 1, it means that the corresponding tuple does not exist. Then we can construct such a mapping to achieve this:

[0066]

[0067] Essentially, the mapping maps m-tuples to a zero-based order based on lexicographic order.

[0068] For example, let (n1,n2,n3) = (2,3,4) and (x1,x2,x3) = (1,0,3), then the length of the bit array is:

[0069] (2+1)(3+1)(4+1)=60

[0070] For example, if int32 is used as the base type of the bit array, two int32s are required, occupying 2*4=8 bytes of memory.

[0071] The index value corresponding to this tuple is

[0072] 1*1+(2+1)*0+((2+1)(3+1))*3=37

[0073] If only one of the above positions is wrong in step (i, j), then only the 37th position of the corresponding bit array is 1, and the rest are 0.

[0074] Therefore, a mapping can be constructed that maps a set of finite edit information to a bit array. Let X be an (N′,m) array, representing the matrix of N′ m-tuples mentioned above. Using vectorized calculation,

[0075]

[0076] is an array of length N′.

[0077] In this way, all possible errors from the first i position of the first string to the first j position of the second string are stored through each bit array. And when the same error occurs at a certain position in the future, all errors can be updated synchronously (the so-called "update" is at the corresponding coordinate +1), and if overflow occurs, it will automatically disappear and will not be updated during the state transfer process.

[0078] Then, based on the mapping described above, the feature mask and translation amount used to update the state are determined. Here, the role of the feature mask is to increase the corresponding coordinates of the comparison codeword by one if an error situation recorded by the feature sequence occurs during the state transition process of dynamic programming, and the generation method is as follows:

[0079] For the mask of the i-th coordinate, a bit sequence is obtained by the above formula using the following set of m-tuples:

[0080] {(x1,…,x m ):x k =0,k≠i,x i =0,1,…,n i -1}

[0081] For the translation amount, it can be obtained according to the following formula

[0082]

[0083] So for a tuple (x1,…,x m ), we can increase the coordinate of the ith position by 1 and keep the rest unchanged.

[0084] For example, let (n1, n2, n3) = (2, 2, 2), the bit array has 27 bits. For each error correction feature, the corresponding error correction mask is: 110110110 11011011 01101101 10 111111000 11111100 01111110 00 111111111 11111111 10000000 00

[0088] The corresponding translation amounts are: 1, (2+1)=3, (2+1)*(2+1)=9.

[0089] Therefore, in the error correction code table generation method according to the embodiment of the present application, determining the feature sequence based on the error correction feature sequence includes: projecting the feature sequence pair corresponding to the unit feature in each error correction feature unit in the error correction feature sequence into the feature sequence.

[0090] Furthermore, in the error correction code table generation method according to the embodiment of the present application, determining the feature mask and the translation amount based on the error correction feature sequence includes: determining an array of an allowable upper limit based on the error correction feature sequence; determining an edit information tuple based on the array of the allowable upper limit; determining a bit array index mapping based on the edit information tuple and the array of the allowable upper limit, the bit array index mapping being used to map the edit information tuple to an order starting from zero based on a lexicographic order; and determining the feature mask and the translation amount based on the bit array index mapping.

[0091] Step S130, iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount.

[0092] Here, during iterative initialization, as described above, the size of the code table is set to n×l, then n empty strings are initialized as code words and stored in the code table, and n(n-1) / 2 dynamic programming matrices are initialized and stored in the matrix array residual, and n(n-1) / 2 0s are initialized to the starting row number of the array.

[0093] Specifically, for any two codewords in the code table, for example, str1 and str2, which represent two strings of equal length, the corresponding dynamic programming matrix is ​​initialized. The initialization process is: first traverse the row and column indexes of the dynamic programming matrix, save the corresponding position with a bit array of length N and all values ​​0, where the size of the dynamic programming matrix is ​​(depthBackTrack+1,2*widthWindow+1); then assign 1 to the 0th position of the bit array of the 0th row and the middle column (widthWindow) of the dynamic programming matrix.

[0094] Here, depthBackTrack and widthWindow represent the backtracking depth and window size respectively. The backtracking depth is that when an error occurs, the character string of the codeword to be generated needs to be matched with the feature sequence pair. If the length of the feature sequence pair is greater than one, the current position of the input character string fragment needs to be backtracked a certain length to match the feature sequence pair. Therefore, the longest feature sequence pair needs to be selected from the error correction feature sequence, and the backtracking number cannot be less than this value.

[0095] In addition, the window size cannot be smaller than this value:

[0096]

[0097] This is because the main factor affecting the sequence synchronization structure is the inconsistency of the amount of change before and after the sequence changes.

[0098] Furthermore, preferably, the window size does not exceed the following value:

[0099]

[0100] This is because this value is the upper limit of the sequence synchronization affected by insertions and deletions.

[0101] That is, in the example where the error types of error correction are limited to include insertion, deletion, and substitution, the number of columns of the dynamic programming matrix is ​​equal to 2 times the maximum number of insertions and deletions plus 1, and the number of rows is equal to the maximum number of callback bits plus 1.

[0102] In one example, the initialization algorithm of the dynamic programming matrix is ​​implemented as follows:

[0103]

[0104] Then, the dynamic programming matrix is ​​iterated based on the ECG algorithm (Edit Computational Graph), where the ECG algorithm is used for a directed acyclic graph (DAG) to determine whether two codewords can be edited from one to another under a given error correction feature sequence.

[0105] Specifically, the iteration of the dynamic programming matrix includes the following steps:

[0106] Step 1: Input two code words in the code table corresponding to the dynamic programming matrix, for example, recorded as seq1 and seq2; the starting coordinates and ending coordinates of the code word comparison, for example, recorded as startIndex and endIndex; the initialized dynamic programming matrix, for example, recorded as mat; the error correction feature sequence, for example, recorded as C; and the starting row number used to mark which row to start the calculation from, for example, recorded as startRow.

[0107] It is worth noting that the algorithm calculates the number of all possible edits from seq1 to seq2 according to C, so the input order of seq1 and seq2 has an impact on the result. However, if C is symmetric, the input order of seq1 and seq2 has no effect.

[0108] Step 2, according to the input matrix size, determine the parameters widthWindow, depthBackTrack; and respectively obtain the feature sequence, for example, recorded as chSeq, the feature mask, for example, recorded as mask and the translation amount, for example, recorded as shiftAmount, and the order of the feature sequence, the feature mask and the translation amount is consistent with the unit order.

[0109] Step 3: Enter the loop. For example, if the loop index is i, then based on the starting coordinate and the ending coordinate, from startIndex to endIndex-1, execute:

[0110] 1) Let d = min(i, widthWindow), that is, the maximum value that needs to be traced back when the string sequence of the first codeword is compared at the i-th position;

[0111] 2) Enter the string sequence of the second code word, index j, the range of j is from -2d to 0:

[0112] I) determining the string sequence indexes index1 and index2 of the first code word seq1 and the second code word seq2 to be compared according to i and j;

[0113] II) Determine whether seq1[index1] == seq2[index2]:

[0114] If yes, use the matrix update algorithm to copy the data of matrix mat about index1-1, index2-1 to index1, index2;

[0115] Otherwise, traverse each feature sequence pair of all error correction feature units in the error correction feature sequence. If seq1[index1], seq2[index2] meet the feature sequence pair, use the matrix update algorithm to update the dynamic programming matrix according to the corresponding mask and translation, otherwise continue to the next feature sequence pair.

[0116] Step 4: Determine the starting row number of the next step startRow = (startRow + endIndex - startIndex) mod (depthBackTrack + 1) according to endIndex and startIndex. Here, the modulus is taken to ensure that startRow does not exceed the number of rows of mat, and return the dynamic programming matrix mat and the starting row number startRow.

[0117] Here, in the loop, index1 and index2 are implemented using the following formula:

[0118] index1=(i+(j / 2))x(j mod 2)

[0119] index2=(i+(j / 2))*((j+1)mod 2)

[0120] In one example, the ECG algorithm is implemented as follows:

[0121]

[0122] It should be noted that lines 1-5 of the above algorithm connect the codewords to be matched in a certain order, and use the indexes i and j to correspond to the indexes of the two codewords. Then, first determine whether the two characters at the corresponding positions are equal, as shown in lines 6-8. If they are equal, directly put the data of the previous line into this row unchanged, otherwise update it according to the feature sequence. In order to be used for GPU computing, the overall algorithm uses bit operations to implement conditional judgments in lines 8, 9 and 13 to avoid thread divergence problems in the GPU.

[0123] Furthermore, the matrix updating algorithm as described above comprises the following steps:

[0124] Step 1: Input the dynamic programming matrix mat, sequence seq1, seq2, corresponding index index1, index2, error correction feature sequence C, error correction mask mask, shift amount shiftAmount, starting row and column rowStart, colStart.

[0125] Here, the starting row and starting column are generally input as 1+startRow–startIndex and widthWindow. This is because the algorithm mainly calculates the offset, so a pair of initial values ​​is needed to clarify the specific rows and columns of the matrix.

[0126] Step 2: According to index1, index2, rowStart, colStart, and the conversion algorithm from index to matrix row and column, get the target row and column rowTarget, colTarget (i.e. the updated matrix coordinates). In addition, in order to prevent the number of rows from overflowing, the number of rows of the matrix needs to be modulo.

[0127] Step 3: Traverse all feature sequence pairs in each error correction feature unit of the error correction feature sequence and enter the loop:

[0128] a) According to the length of the two feature strings, calculate the new indexes index1New and index2New respectively. The above indexes will be used as the starting indexes of the fragments to be matched. At the same time, calculate the corresponding source row and column numbers rowSource and colSource;

[0129] b) Determine whether the error occurs:

[0130] If it matches, perform a bitwise AND operation on the data of mat[rowSource,colSource] and mask, shift shiftAmount bits to the higher bit, and perform a bitwise AND operation on mat[rowTarget,colTarget]. Specifically, the bitwise AND operation is equivalent to taking a value, the shift is equivalent to updating the source data, and the bitwise OR operation is equivalent to adding the result to the target position;

[0131] If not satisfied, the next feature is matched.

[0132] Here, the matrix update algorithm only works on the dynamic programming matrix mat, so there is no need to return the result.

[0133] In one example, the matrix update algorithm is implemented as follows:

[0134]

[0135]

[0136] It should be noted that in the above matrix update, the dynamic programming matrix mat is updated by performing an AND operation on the corresponding mask and the corresponding position data of the starting row and column. This process is code acquisition. Then after high-bit translation, the data corresponding to the position of the target row and column is assigned through an OR operation. Therefore, for a given error correction feature sequence, the above process is equivalent to adding +1 to the specific coordinates of all edit records, and the other coordinates remain unchanged. If you want to copy the data at a certain position directly, then the mask must be all 1 data, and the shiftAmount must be 0. If you choose not to perform the update of the error correction feature sequence, you can set the mask to 0, so no matter how it is translated, the data of the target row and column will remain unchanged.

[0137] Furthermore, the algorithm for converting the above indicators to matrix rows and columns is as follows:

[0138] Step 1: Input index1, index2, rowStart, colStart;

[0139] Step 2: Get the column offset bias = index2-index1;

[0140] Step 3, calculate the offset of the row i = (index1 + index2 + abs(bias)) / 2; and

[0141] Step 4, calculate the row and column: row = i + rowStart, col = bias + colStart.

[0142] In one example, the algorithm for converting the above indexes to matrix rows and columns is implemented as follows:

[0143]

[0144] In addition, judging whether the sequence segment satisfies a certain characteristic sequence pair in the error correction characteristic unit includes:

[0145] Step 1: Input segments seq1, seq2, indexes index1, index2 to be matched, and feature sequence c = (p1, p2);

[0146] Step 2: First determine whether the length of seq1 is equal to the length of p1, and do the same operation on seq2 and p2;

[0147] Step 3: Traverse the characters of seq1 and enter the loop:

[0148] a) Determine whether the corresponding position of p1 is '*'. If so, it means it meets the requirements;

[0149] b) Consider whether there is a '=':

[0150] I) If yes, check whether the characters in seq1 are equal to the characters in the corresponding positions of seq2. If they are equal, then the match is met; otherwise, the match is not met;

[0151] II) If not, consider whether the character in seq1 is equal to the character in the corresponding position of p1. If they are equal, then it is in compliance; otherwise, it is not in compliance;

[0152] Step 4: Repeat the above loop and judgment for seq2 and p2;

[0153] Step 5: Finally, when all judgments are met, return true (True); otherwise, return false (False).

[0154] In one example, the algorithm for determining whether the error at the current sequence position satisfies the condition of the error correction feature unit is implemented as follows:

[0155]

[0156] That is, first determine whether the lengths of the two strings to be compared and the corresponding feature strings are consistent; then determine whether the two feature strings contain wildcards; if not, determine whether there is an equal sign. If yes, determine whether the corresponding position of the first string to be compared is equal to the corresponding position of the second string to be compared; otherwise, determine whether the corresponding positions of the first string to be compared and the first feature string are equal, and if they are not equal, return zero. The same is true for the second string to be compared. If all of the above are satisfied, return one.

[0157] At this time, the update of the dynamic programming matrix is: for the index1, index2 of the sequence seq1, seq2, as well as the current feature sequence pair (p1, p2), mask mask and translation amount, use mask to take out the data of the matrix at the index coordinate (index1-|p1|, index2-|p2|), shift the amount to the high position, and then put it at the position of the matrix (index1, index2) through the or operation. Therefore, this process is equivalent to adding +1 to a specific coordinate by bit operation. This is because a group of edit count tuples are bit-encoded, and each 1 position in the bit array corresponds to a different count tuple. When a certain position of all count tuples in the group is +1, it is equivalent to a certain amount of translation on the bit array.

[0158] That is, in the error correction code table generation method according to the embodiment of the present application, based on the coding principle as described above, the dynamic programming matrix is ​​updated by the incremental calculation method of "bottom-up", that is, the dynamic programming matrix of the current step is retained for each calculation, and the next calculation is performed directly on the result of the matrix. In addition, by limiting the matching range of the codeword, only the error types that need to be concerned and the corresponding number of errors are determined, that is, the number of errors that may occur in two codewords is limited, so only the string matching of the codeword within a specific range needs to be considered. For example: if the maximum number of insertions and deletions that need to be concerned is 2, then the strings separated by more than 2 bits do not need to consider whether they match, so computing power is saved. In addition, for the string fragments that need to be investigated, their length will affect the backtracking depth that needs to be confirmed forward for each match.

[0159] Among them, the size of the dynamic programming matrix is ​​(maximum number of callback bits + 1, 2 times the maximum number of insertions and deletions + 1), and the (i mod number of rows, j) of the matrix stores the remaining editing quota consumed by codeword s1[:i] to reach codeword s2[:j] under a given "error scenario". Therefore, starting from the (0,0) position, for any (i,j) position, first determine whether s1[i] is equal to s2[j]. If they are equal, it means a match, so the data of (i-1,j-1) is directly copied to the position of (i,j). If they are not equal, the state transfer is performed according to the matrix update method described above, and the data is updated to the position of (i,j). Note that in the one-to-one correspondence algorithm between the row and column numbers of the above matrix and the index1 and index2 of the sequence, only the row mod the number of rows is used. That is, the bit array obtained according to the error correction feature sequence is bitwise ANDed under the corresponding matrix element according to the callback number, and then the result is bitwise shifted to the high position, the overflow part is discarded, the low position is filled with zero, and finally the result is bitwise ORed and accumulated to the current array. When performing this step, it is necessary to determine whether the position where the error occurs in the current s1 sequence meets the error correction feature. If it does, continue, otherwise skip. If all positions are calculated, the algorithm stops.

[0160] Therefore, in the error correction code table generation method according to the embodiment of the present application, iterating the dynamic programming matrix based on the feature sequence, the feature mask and the translation amount includes: when the size of the code table is n×l, initializing n empty strings as the code words of the code table, and the length of each code word is l; initializing n(n-1) / 2 dynamic programming matrices, wherein each dynamic programming matrix corresponds to two code words of the code table, and its size is (backtracking depth + 1, 2*window size + 1), each position of the dynamic programming matrix is ​​a bit array, and the 0th position of the bit array of the 0th row middle column of the initialized dynamic programming matrix is ​​1; and, initializing n(n-1) / 2 starting rows.

[0161] Furthermore, in the above-mentioned error correction code table generation method, iterating the dynamic programming matrix based on the feature sequence, the feature mask and the translation includes: iterating the dynamic programming matrix using an edit computation graph algorithm, wherein in the i-th iteration, the starting row number is determined for the starting coordinates and the ending coordinates of the codeword comparison, the comparison index of the two codewords is determined based on i and the backtracking depth, and the comparison string pair of the two codewords based on the comparison index is determined; in response to the two comparison string pairs being equal, updating the dynamic programming matrix using a first matrix update algorithm; and in response to the two comparison string pairs being unequal and conforming to the feature sequence, updating the dynamic programming matrix using a second matrix update algorithm according to the feature mask and the translation.

[0162] In addition, in the above-mentioned error correction code table generation method, using the first matrix update algorithm to update the dynamic programming matrix includes: determining the coordinates of the dynamic programming matrix with the comparison index of the two code words minus one respectively; taking out the data of the dynamic programming matrix at the coordinates with a mask; shifting one to the high bit and placing it in the position of the dynamic programming matrix with the comparison index of the two code words as the coordinates through an OR operation.

[0163] In addition, in the above-mentioned error correction code table generation method, updating the dynamic programming matrix using a second matrix update algorithm according to the feature mask and the translation amount includes: determining the coordinates of the dynamic programming matrix with the comparison index of the two code words minus the current feature sequence; extracting the data of the dynamic programming matrix at the coordinates with a mask; shifting the amount of the translation to a high position; and, placing it in the position of the dynamic programming matrix with the comparison index of the two code words as the coordinates through an OR operation.

[0164] Step S140, in response to the loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, determining two code words in a code table corresponding to each dynamic programming matrix completed by iteration to obtain an error correction code table.

[0165] For the array of bit arrays in the dynamic programming matrix, the loss function is defined as:

[0166] min{i:array[i]==1,0≤i <len(array)}

[0167] Then for the entire dynamic programming matrix mat, its loss function is equal to

[0168]

[0169] Wherein d is the backtracking depth, w is the window size, s is the starting function, and (i, j) is the comparison index of the two codewords.

[0170] Here, because the loss function of each element in the matrix is ​​equivalent to the upper bound of the maximum edit count of the corresponding parts of the estimated codewords seq1 and seq2, the result of summing and averaging each element in the matrix can reflect how many ways the two codewords can be changed from one to the other (of course, this edit number is constrained according to the error correction feature sequence). Therefore, if the loss function is larger, it means that there are many ways to edit, that is, the more similar; and the smaller the loss function, the less similar. Therefore, in the error correction code table generation method according to the embodiment of the present application, the loss function is less than the predetermined threshold value, indicating that the task is completed.

[0171] That is, if the loss function of the dynamic programming matrix is ​​less than or equal to a predetermined threshold, for example, the value of the loss function is zero, then it means that the distinction between the two code words seq1 and seq2 in the array code table corresponding to the dynamic programming matrix is ​​sufficient to achieve the error correction goal, and no confusion will be caused by adding a sequence backward.

[0172] In this way, if n(n-1) / 2 dynamic programming matrices complete the iteration, an error correction code table containing n code words can be obtained.

[0173] Therefore, in the error correction code table generation method according to the embodiment of the present application, the loss function of the dynamic programming matrix is ​​expressed as:

[0174]

[0175] Wherein d is the backtracking depth, w is the window size, s is the starting function, and (i, j) is the comparison index of the two codewords.

[0176] Moreover, the ECG algorithm can be degenerated into the calculation of the similarity of the previous sequence. Specifically, in the error correction feature sequence, a sufficient upper limit of the number of errors is set, and the feature sequence can be set to the corresponding editing mode. According to the above indicator to matrix row and column conversion algorithm, let index1 be the length of s1, index2 be the length of s2, so as to obtain the relative column coordinate col, then according to the final startRow, the loss function value Loss(mat[startRow][col]) of the dynamic programming matrix mat at startRow, col can be obtained, which is the corresponding editing distance.

[0177] Here, for the iteratively completed dynamic programming matrix, the bit array at the position of Loss = mat[startRow][widthWindow] has the following information:

[0178] If the bit array is all 0, it means that the two codewords to the current state have sufficient discrimination (the discrimination is defined by the error correction feature sequence);

[0179] The number of 1s in the bit array reflects the number of all possible editing methods from the first codeword to the second codeword under consideration;

[0180] The size of the index value with the lowest bit being 1 can be used to reflect the extent to which two code words have sufficient discrimination, that is, the bit array is arranged from low to high, and the index value with the lowest bit being 1 is recorded. The smaller the index value, the lower the discrimination.

[0181] min{i:Loss[i]==1}

[0182] That is, the bit array is initialized according to the error correction feature sequence as described above, where the number of bits of the bit array is the product of the maximum number of errors, and different consecutive blocks represent error situations with different error correction features. In the iteration of the dynamic programming matrix, if the "error quota" is not consumed, the default lowest bit is 1 and the remaining bits are 0; if it is consumed, all bits are 0. 1 in different positions represents the type of error that may occur. Therefore, the loss function is essentially defined as the distance from the lowest bit 1 to the highest bit in the bit array of the row and middle column where the dynamic programming matrix stays at the end of the algorithm. (For example: 0010, loss = 3; 0000, loss = 0). Therefore, when the loss functions of all two codeword sequences in the code table are 0, it means that the code table has the error correction capability to complete the corresponding error correction task.

[0183] Here, when the setting of the error correction feature sequence conforms to the traditional editing method and the error upper limit is set high enough, the error upper limit value will be equal to the corresponding editing distance, including Hamming distance, LCS distance, Levenshtein distance, Damerau-Levenshtein distance, etc. Among them, each of the above editing distances corresponds to different error correction requirements, for example, the Hamming distance only considers the minimum number of editing for substitution; the LCS distance only considers the minimum number of editing for the sum of insertion and deletion; the Levenshtein distance considers the minimum number of editing for the sum of substitution, insertion and deletion; and the Damerau-Levenshtein distance considers the minimum number of editing for the sum of insertion, deletion, substitution and adjacent swapping.

[0184] In addition, in the error correction code table generation method according to the embodiment of the present application, incremental calculation can be used to iterate the dynamic programming matrix to generate the code table. Specifically, when initializing the dynamic programming matrix mat, widthWindow=mat.shape[0] / / 2, depthBackTrack=mat.shape[1].

[0185] Therefore, if there are sequences seq1+seg1 and seq2+seg2 that need to be incrementally processed, it can be determined that the two input sequences are

[0186] seq1'=(seq1+seg1)[-widthWindow-|seg1|:],

[0187] seq2'=(seq2+seg2)[-widthWindow-|seg1|:].

[0188] Here, by default, |seg1|=|seg2|. If they are not equal, spaces are added to the shorter segment to serve as placeholders. Note: s[-a:] means all characters from max{0,len(s)-a} to the end.

[0189] At the same time, startIndex and endIndex can be determined as follows:

[0190] endIndex=len(seq1'), startIndex=endIndex-len(seg1).

[0191] Therefore, the incremental algorithm can be determined as follows: Compared with the above ECG algorithm, the incremental calculation determines startIndex and endIndex according to the size of the newly added fragment, and then further calculates the previously calculated dynamic programming matrix mat, that is, there is no need to iterate the dynamic programming matrix from the beginning until the new dynamic programming matrix mat and the starting number of rows startRow are returned.

[0192] In one example, the ECG delta calculation algorithm is implemented as follows:

[0193]

[0194] Therefore, in the error correction code table generation method according to the embodiment of the present application, it further includes: when the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, the starting coordinates and the ending coordinates of the dynamic programming matrix are determined based on the length of each codeword in the new code table and the window width during iteration; and, based on the starting coordinates and the ending coordinates, a new starting number of rows is determined, and based on the new starting number of rows, the dynamic programming matrix is ​​continued to be iterated on the basis of the dynamic programming matrix corresponding to the generated code table to obtain the new code table.

[0195] Schematic encoding and decoding method

[0196] Figure 2 FIG. 1 is a flowchart of a coding method based on a generated code table according to an embodiment of the present application. Figure 2 As shown, the encoding method based on generating a code table according to an embodiment of the present application includes the following steps.

[0197] Step S210, determining the floor value of the logarithmic function with base 2 of the code table size of the generated code table. That is, assuming that the code table size of the generated code table is n, determining

[0198] Step S220, pre-processing the input binary sequence to pad the end with zeros and recording the number of zeros padded.

[0199] Step S230, group encoding the pre-processed input binary sequence with the number of the rounded-down values ​​as a group. That is, assuming that the length of the input binary sequence is N, The bits are grouped into a group for group coding, and the code corresponding to each n′-length binary sequence is the decimal number corresponding to the binary sequence in the code table.

[0200] That is, because the group encoding is performed based on the number of the rounded-down values, the input binary sequence needs to be preprocessed to pad the end with zeros and record the number of zero padding. Here, for the input binary sequence that has completed the preprocessing, the number of zero padding will not exceed n′, so it is necessary The number of zeros to be padded is recorded in binary format. The number of zeros to be padded is as shown in the following formula:

[0201] s=(n′-(m+N)mod n′)mod n′.

[0202] In this way, for the preprocessed binary sequence, its length is divisible by n', so it can be divided into equal parts, each of which is n' long, and encoded according to the corresponding method of the code table, and finally the encoded string sequence is spliced ​​together in order to finally obtain a string segment. In addition, the length of the obtained string segment should be an integer multiple of the length L of the code word in the code table.

[0203] Figure 3 FIG. 1 is a flowchart of a decoding method based on a generated code table according to an embodiment of the present application. Figure 3 As shown, the decoding method based on the generated code table according to the embodiment of the present application includes the following steps.

[0204] Step S310, the received sequence to be decoded is divided into a predetermined number of groups, the predetermined number is the number of coding groups, and the length of each group is the length of each codeword of the generated code table. That is, for the received sequence to be decoded, assuming that the code table size of the generated code table is n, n'=nlog2n is determined, and the sequence to be decoded is divided into n' groups, and the length of each group is the length L of the codeword of the generated code table.

[0205] Step S320, starting from the first group, each group is compared with the code table to determine the string corresponding to the group in the code table. That is, starting from the first group seq[0:L], it is compared with each element in the code table. If the current error does not exceed the range, there is only one sequence corresponding to it, that is, the sequence is correct, and the next step is entered; otherwise, the over-range processing will be performed.

[0206] Step S330, for the kth group in the sequence to be decoded and the corresponding error correction sequence index is i, the codeword in the code table is added with the prefix of the codeword minus the traceback depth, and matched with the sequence of the previous correction sequence minus the start of the traceback depth and the length of the codeword multiplied by k.

[0207] Specifically, for the kth group in the sequence to be decoded and the corresponding error correction sequence index is i, the prefix codebook[i][-widthWindow:] is added to the sequence in the code table, that is, the i-th codeword in the code table minus the backtracking depth, and it is matched with the string of the previous sequence minus the backtracking depth, that is, seq[L*(k-1)–widthWindow:L*k]. That is, considering the kth segment (k>1) and the case where the last corrected sequence index is i, the prefix codebook[i][-widthWindow:] is added to the sequence in the code table, and it is matched with seq[L*(k-1)–widthWindow:L*k]. If the current range is not exceeded, there is only one correct sequence corresponding to it, and step S340 is entered; otherwise, step S350 is entered.

[0208] Step S340: After all groups of the sequence to be decoded are matched, all correct sequences are concatenated to obtain a decoded sequence.

[0209] In addition, in the decoding method based on the generated code table according to the embodiment of the present application, an over-range processing is further included, including: when the string corresponding to the group is all 0, the last sequence of zero is taken to correspond to it, and the subsequent default is the first codeword codebook[0] in the code table; when there are multiple non-0 strings corresponding to the group, the sequence with the largest loss function is taken to correspond to it. Here, the loss function is the above loss function:

[0210]

[0211] That is, because the loss function of each element in the matrix is ​​equivalent to estimating the upper bound of the maximum edit count of the corresponding part of codeword s1 to codeword s2, the result of summing and averaging can reflect how many ways there are for two codewords to change from one to the other (of course, this edit number must be constrained according to the error correction feature sequence). If the loss function is larger, it means there are many ways to edit, that is, the more similar; and the smaller the loss, the less similar. In the algorithm for code table generation, the smaller the Loss, the more it means that the task is completed. But on the contrary, in decoding, contrary to what is described above, the larger the loss function, the more it means that it is the correct sequence.

[0212] Schematic method for generating gene sequence code table

[0213] Figure 4 The figure illustrates a schematic flow chart of a method for generating an error correction code table for a gene sequence according to an embodiment of the present application.

[0214] like Figure 4 As shown, the method for generating an error correction code table for a gene sequence according to an embodiment of the present application includes the steps of:

[0215] S410, determining the size of the error correction code table for the gene sequence, and generating an error correction feature sequence based on the gene sequence {A, T, C, G} as a finite character table and the error correction error types and error correction error numbers corresponding to the gene sequence error correction target.

[0216] That is, in the storage of gene sequences, such as DNA, the error correction code table generation method described above can be used to design an error correction code table for the gene sequence, because the gene sequence space can also be regarded as a data coding space, and the received gene sequence is evaluated by a suitable function or algorithm to determine whether there is an error and attempt to correct it.

[0217] In addition, there are two major challenges in existing DNA storage coding schemes: different errors occur simultaneously and unevenly; and the structure of the sequence affects decoding.

[0218] First, the imbalance of errors is manifested in that in the process of artificially synthesizing DNA, there are many types of errors (such as insertion, deletion, replacement, simultaneous continuous insertion and deletion, swap, replication, etc.) and they will occur at the same time. At the same time, for different bases, the frequency of occurrence of different error modes is also different. Despite this, they can all be unified as insertion + deletion errors, and then error correction can be implemented based on VT codes (most of the work is done this way). However, this often introduces more redundancy. For example, for 1 substitution error, the error correction code is often required to have the error correction capability of 1 insertion + 1 deletion; for 1 swap, the error correction capability of 2 insertions + 2 deletions is required. To this end, error correction design is required for different bases and different error types, but it is quite difficult to achieve this by designing rules.

[0219] In addition, the stability of the sequence structure needs to be taken into account. In other words, in DNA storage, another difference from traditional encoding is that the structure of the sequence needs to be considered. From synthesis to sequencing, if the GC content of the DNA sequence deviates too much, and the frequency of repeated occurrence of the same base is too high, it will seriously affect the final result. The balance of GC content is generally called GC-Balance; the maximum frequency of repeated occurrence of the same base is called runlength.

[0220] Therefore, using the error correction code table generation method as described above, when determining the length of each group of gene sequences to be encoded, the different permutations and combinations of 01 are determined, so that, through the encoding strategy based on the above error correction code table, each binary sequence has a unique base sequence corresponding to it. Then, such a base sequence code table can be constructed to achieve the corresponding error correction requirements, and once the code table is obtained, it can be encoded by the index of the encoding position, and decoding and error correction can be achieved through comparison.

[0221] Step S420, determining a feature sequence, a feature mask and a shift amount based on the error correction feature sequence.

[0222] Step S430, iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount.

[0223] Step S440, in response to the loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, determining two code words in a code table corresponding to each dynamic programming matrix completed by iteration to obtain an error correction code table for a gene sequence.

[0224] In the above-mentioned error correction code table generation method for gene sequences, the error correction feature sequence is generated based on the gene sequence {A, T, C, G} as a finite character table and the error correction error type and the number of error correction errors corresponding to the gene sequence error correction target, including: generating a complete character table based on the gene sequence {A, T, C, G} as a finite character table; defining a feature string space by a non-intersecting union of multiple complete character tables, the feature string space containing an empty character string representation; defining a feature sequence pair representing the error correction error type by two different elements in a monoid defined in the feature string space by character splicing; defining an error correction feature unit based on multiple feature sequence pairs and an allowable upper limit that is a positive integer greater than or equal to zero; and defining an error correction feature sequence based on a sequence composed of multiple error correction feature units.

[0225] In the above-mentioned method for generating an error correction code table for a gene sequence, the characteristic sequence pairs representing the types of error correction errors include:

[0226] represents the replacement (p1, p2), p i ≠-,i=1,2;

[0227] represents the insertion of (-,p2), p2≠-; and

[0228] Represents deleted (p1,-), p1≠-;

[0229] Where p1, p2 are two different elements in the monoid, and - is the empty string.

[0230] In the above-mentioned error correction code table generation method for gene sequences, after generating the error correction feature sequence, it further includes symmetrizing the error correction feature sequence, including: starting to traverse each error correction feature unit for the error correction feature sequence; multiplying the allowable upper limit of each error correction feature unit by two; and, for each feature sequence pair in the error correction feature unit, determining whether a symmetric feature sequence pair that is arranged opposite to two elements of the feature sequence is in the error correction feature unit, and in response to the symmetric feature sequence pair not being in the error correction feature unit, placing the symmetric feature sequence pair in the error correction feature unit.

[0231] In the above-mentioned error correction code table generation method for gene sequences, determining the feature sequence based on the error correction feature sequence includes: projecting the feature sequence pair corresponding to the unit feature in each error correction feature unit in the error correction feature sequence into the feature sequence.

[0232] In the above-mentioned error correction code table generation method for gene sequences, determining the feature mask and the translation amount based on the error correction feature sequence includes: determining an array of an allowable upper limit based on the error correction feature sequence; determining an edit information tuple based on the array of the allowable upper limit; determining a bit array index mapping based on the edit information tuple and the array of the allowable upper limit, the bit array index mapping being used to map the edit information tuple to an order starting from zero based on a lexicographic order; and determining the feature mask and the translation amount based on the bit array index mapping.

[0233] In the above-mentioned error correction code table generation method for gene sequences, iterating the dynamic programming matrix based on the feature sequence, the feature mask and the translation includes: when the size of the code table is n×l, initializing n empty strings as the code words of the code table, and the length of each code word is l; initializing n(n-1) / 2 dynamic programming matrices, wherein each dynamic programming matrix corresponds to two code words of the code table, and its size is (backtracking depth+1, 2*window size+1), each position of the dynamic programming matrix is ​​a bit array, and the 0th position of the bit array of the 0th row middle column of the initialized dynamic programming matrix is ​​1; and, initializing n(n-1) / 2 starting rows.

[0234] In the above-mentioned error correction code table generation method for gene sequences, iterating the dynamic programming matrix based on the feature sequence, the feature mask and the translation includes: iterating the dynamic programming matrix using an edit computation graph algorithm, wherein in the i-th iteration: for the starting row number determined by the starting coordinates and the ending coordinates of the codeword comparison, the comparison index of the two codewords is determined based on i and the backtracking depth, and the comparison string pair of the two codewords based on the comparison index is determined; in response to the two comparison string pairs being equal, the dynamic programming matrix is ​​updated using a first matrix update algorithm; and, in response to the two comparison string pairs being unequal and conforming to the feature sequence pair, the dynamic programming matrix is ​​updated using a second matrix update algorithm according to the feature mask and the translation.

[0235] In the above-mentioned error correction code table generation method for gene sequences, using the first matrix update algorithm to update the dynamic programming matrix includes: determining the coordinates of the dynamic programming matrix with the comparison index of the two code words minus one respectively; taking out the data of the dynamic programming matrix at the coordinates with a mask; and shifting one to the high position and placing it in the position of the dynamic programming matrix with the comparison index of the two code words as the coordinate through an OR operation.

[0236] In the above-mentioned error correction code table generation method for gene sequences, updating the dynamic programming matrix using a second matrix update algorithm according to a feature mask and a translation amount includes: determining the coordinates of the dynamic programming matrix with the comparison index of the two code words minus the current feature sequence pair; extracting the data of the dynamic programming matrix at the coordinates with a mask; shifting the amount of the translation amount to a high position; and, placing it at a position of the dynamic programming matrix with the comparison index of the two code words as the coordinates through an OR operation.

[0237] In the above-mentioned error correction code table generation method for gene sequences, the loss function of the dynamic programming matrix is ​​expressed as:

[0238]

[0239] Wherein d is the backtracking depth, w is the window size, s is the starting function, and (i, j) is the comparison index of the two codewords.

[0240] In the above-mentioned error correction code table generation method for gene sequences, it further includes: when the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, the starting coordinates and the ending coordinates of the dynamic programming matrix iteration are determined based on the length of each codeword in the new code table and the window width; and, based on the starting coordinates and the ending coordinates, a new starting row number is determined, and based on the new starting row number, the dynamic programming matrix is ​​continued to be iterated on the basis of the dynamic programming matrix corresponding to the generated code table to obtain the new code table.

[0241] Here, those skilled in the art can understand that the specific details of the above error correction code table generation method for gene sequences have been described in detail in the previous "Schematic code table generation method" section, and will not be repeated here to avoid redundancy.

[0242] As mentioned above, for encoding in DNA storage, it is necessary to consider the impact of the structure of the sequence on the stability of the actual DNA sequence, such as the maximum number of consecutive repetitions of the same base and the proportion of the sum of the GC content to the whole. Therefore, constraints should be set for each increment. That is, for the incremental strategy, the probability of selecting a base for each sequence in the next step is actually returned by inputting the information of the current code table, and then randomly selecting according to the probability. Therefore, the sequence structure constraint is further adjusted by a set of probabilities returned by the incremental strategy, and the process must also satisfy the requirement of being differentiable.

[0243] Here, each input is a code table and the current state of the code table. Assume that the code table includes n code words, and the length of each code word is l, and the size of the finite alphabet is m (usually m=4). Each base can have the following correspondence:

[0244] A~(1,0,0,0)=e1

[0245] G~(0,1,0,0)=e2

[0246] C~(0,0,1,0)=e3

[0247] T~(0,0,0,1)=e4

[0248] Therefore, the code table is an array of size [n, l, 4], and the code table state should be of size [n, 4+1+1], where for the code table, the first indicator represents the id of the code word, the second indicator represents the id of the sequence base, and the third indicator is used to correspond to a specific base; for the code table state, the first indicator represents the sequence id, the first 4 digits of the second indicator are used to record consecutively appearing characters, the next digit is used to record the number of consecutive appearances, and the last digit records the total GC content.

[0249] Assume that the probability of each base increment in a sequence at the current step is (x1, x2, x3, x4), the current number of repetitions is recorded as run, and the current total GC content is recorded as g. Then it is agreed that the maximum number of consecutive identical bases does not exceed L r The difference between GC ratio and 50% is gc var .

[0250] The operators for GC ratio constraints are:

[0251] AT prob =x1+x4,GC prob=x2+x3

[0252] AT′ prob =(AT prob *ReLU(g-(0.5-gc var ))+(1-2AT prob )*ReLU(g-0.5)+(AT prob -1)*ReLU(g-(0.5+gc var ))) / gc var

[0253] GC′ prob =1-AT′ prob

[0254] The probability of a new base is:

[0255]

[0256] Where ReLU(x)=max{x,0}.

[0257] This operator can ensure that when GC exceeds 0.5+gc var When x2 and x3 are 0, the value is 0; otherwise, if it is lower than 0.5-gc var When , the values ​​of x1 and x4 are 0. When it is equal to 0.5, the original distribution ratio is not changed.

[0258] The operator for consecutive identical base constraints is:

[0259] After GC constraint, the obtained AGCT probability is recorded as (x1, x2, x3, x4).

[0260]

[0261] This ensures that as run increases, the probability of consecutive bases appearing again decreases; when run equals the upper limit, the probability drops to 0.

[0262] Therefore, in the method for generating an error correction code table for a gene sequence according to an embodiment of the present application, before determining the starting coordinates and the ending coordinates of the dynamic programming matrix iteration based on the length of each codeword and the window width in the new code table, the method further includes:

[0263] When the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, a base sequence is added as a suffix to each generated codeword according to a predetermined probability corresponding to each gene base based on an incremental strategy.

[0264] Furthermore, in the above-mentioned error correction code table generation method for gene sequences, the increment strategy is determined based on the suffix increment constraint corresponding to the gene sequence, and the suffix increment constraint includes a GC ratio constraint and a continuous identical base constraint.

[0265] In addition, the incremental strategy can also be determined based on other strategies, such as the Monte Carlo strategy, that is, when adding suffixes in each round, the probability of each base is 0.25, and the probability is incrementally constrained according to the previously generated sequence. For example, 1-10 bases are added continuously (predetermined) and repeated for 10-15 rounds, and the set of results with the largest reduction in the loss function is selected and returned.

[0266] In addition, the incremental strategy can also be determined based on deep reinforcement learning. For example, the code table is first processed as follows:

[0267]

[0268] Therefore, X is an array of size [N,4].

[0269] Input it into the deep neural network and return an array Y of size [N,4], which means the probability of adding different bases to each sequence.

[0270] The suffix enhancement constraint, such as the GC constraint and the operator of the continuous same constraint mentioned above, will be returned through the state of the current code table, and used to act on Y to finally obtain the result Y′, which is the final probability, and the size is also [N, 4].

[0271] According to the probability Y′, a set of candidate suffixes is generated, and the optimal suffix (with the smallest loss function) is calculated according to the incremental ECG algorithm and added to the end of the code table. At the same time, the current reward and action are recorded.

[0272] As the loss function decreases until it reaches 0, the algorithm stops and the policy is optimized through policy gradients based on the rewards and actions obtained.

[0273] In one example, the rewards are designed as follows:

[0274] First, use the Monte Carlo method to obtain a code table, record its length L0, and record R0 of the loss function.

[0275] Then we enter the reinforcement learning phase. According to the suffix added at each i-th step, there will be a new loss function l, which is compared with the loss function of the corresponding step of R0.

[0276] If l≤L0, the current reward is (R0[i]-l) / (R0[i]).

[0277] If l>L0, the current reward is -1.

[0278] If the sequence length in the code table of the trial result is less than L0, record the new length L1 and the record R0 of the loss function. Repeat the above steps.

[0279] That is, in the above-mentioned error correction code table generation method for gene sequences, the incremental strategy is based on the Monte Carlo strategy. In the Monte Carlo strategy, the probability of each base when adding a suffix in each round is 0.25, and the probability is incrementally constrained based on the loss function according to the code words of the generated code table.

[0280] Furthermore, in the above-mentioned error correction code table generation method for gene sequences, the incremental strategy is based on a deep reinforcement learning strategy, and the deep reinforcement learning strategy includes:

[0281] Inputting the code table into a deep neural network to obtain a feature matrix, wherein the feature matrix represents the probability of adding different bases to each code word of the code table;

[0282] Obtaining a result matrix from the feature matrix based on a suffix increment constraint;

[0283] Generate candidate suffixes based on the result matrix, and calculate the optimal suffix, reward and action by minimizing the loss function;

[0284] The deep reinforcement learning strategy is optimized through policy gradient according to the reward and action.

[0285] In addition, in the above-mentioned error correction code table generation method for gene sequences, the reward is:

[0286] Obtain the code table, code table length L0, and loss function R0 through the Monte Carlo method;

[0287] In the reinforcement learning process, the new loss function l after adding the suffix in the i-th step is compared with the loss function of the corresponding step of R0 to determine the reward:

[0288] In response to l≤L0, the reward is (R0[i]-l) / (R0[i]);

[0289] In response to l>L0, the reward is -1.

[0290] Schematic gene sequence encoding and decoding method

[0291] Figure 5 The figure illustrates a schematic flow chart of an encoding method based on a gene sequence code table according to an embodiment of the present application.

[0292] like Figure 5As shown, the encoding method based on the gene sequence code table according to the embodiment of the present application includes the following steps:

[0293] S510, determine the floor value of the logarithmic function with base 2 of the code table size of the gene sequence code table; S520, preprocess the input binary sequence to pad the end with zeros and record the number of padding zeros; and, S530, group encode the preprocessed input binary sequence according to the number of the floor values ​​as a group.

[0294] Here, in the process of storing the gene sequence, since the order of each fragment is shuffled, a part needs to be reserved to retain the sequence order information. Assuming that the size of the code table for generating the code table is n×L, and the length of the input binary sequence is N, the input binary needs to be sorted according to The encoding of each n′ long binary sequence is the decimal number corresponding to the binary sequence index in the code table. Before transcoding, the input binary needs to be preprocessed, that is, the number of zero padding at the end and the data header recording the number of zero padding.

[0295] For the data header, the number of zero padding will not exceed n', so it is necessary The number of zeros to be padded is recorded in binary format, and the number of zeros to be padded is as shown in the following formula:

[0296] s=(n′-(m+N)mod n′)mod n′.

[0297] For the preprocessed binary sequence, its length is divisible by n′, so it can be divided into equal parts, each with a length of n′, and encoded according to the corresponding method of the code table. Finally, the encoded gene fragments are spliced ​​together in order to obtain the gene sequence.

[0298] Therefore, the length of the obtained gene sequence should be an integer multiple of L. Since the length of each gene sequence fragment in the gene sequence storage is constrained (between 100 and 200 bases), there are also certain requirements for the input binary sequence, that is,

[0299]

[0300] Therefore, in the encoding method based on the gene sequence code table according to the embodiment of the present application, when the code table size of the generated code table is n×l and the length of the binary sequence to be encoded is N, N satisfies:

[0301]

[0302] in And s=(n′-(m+N)mod n′)mod n′.

[0303] Figure 6 The figure illustrates a schematic flow chart of a decoding method based on a gene sequence code table according to an embodiment of the present application.

[0304] like Figure 6 As shown, the decoding method based on the gene sequence code table according to the embodiment of the present application includes the following steps:

[0305] S610, dividing the received gene sequence to be decoded into a predetermined number of groups, the predetermined number being the number of coding groups of the logarithmic function with base 2 of the code table size of the gene sequence code table, and the length of each group being the length of each codeword in the gene sequence code table; S620, comparing each group with the gene sequence code table starting from the first group to determine the codeword corresponding to the group in the gene sequence code table; S630, for the kth group in the gene sequence to be decoded and the corresponding error correction sequence index i, adding a prefix of the codeword minus the window size to the codeword in the gene sequence code table, and matching it with a sequence of the length of the previous correction sequence minus the start of the window size and the length of the codeword multiplied by k; and, S640, after matching all groups of the gene sequence to be decoded, splicing all correct sequences to obtain a decoded gene sequence.

[0306] That is, for the kth group in the gene sequence to be decoded and the corresponding error correction sequence index is i, the prefix codebook[i][-widthWindow:] is added to the sequence in the code table, that is, the i-th codeword in the code table minus the backtracking depth, and it is matched with the string of the previous sequence minus the backtracking depth, that is, seq[L*(k-1)–widthWindow:L*k]. That is, considering the kth segment (k>1), and the case where the last corrected sequence index is i, the prefix codebook[i][-widthWindow:] is added to the sequence in the code table, and it is matched with seq[L*(k-1)–widthWindow:L*k].

[0307] In addition, the above decoding method based on generating a code table further includes over-range processing: when the character strings corresponding to the group are all 0, the last sequence of zero is taken to correspond to it, and subsequently defaulted to the first codeword codebook[0] in the code table; and, when there are multiple non-zero characters in the character string corresponding to the group, the sequence with the largest loss function is taken to correspond to it.

[0308] Application Examples

[0309] Figure 7A schematic diagram of an application example of a method for generating an error correction code table for a gene sequence and a coding and decoding method based on a gene sequence code table according to an embodiment of the present application is illustrated.

[0310] like Figure 7 As shown, in the process of generating the code table, the loss function is minimized by adding base suffixes in a "bottom-up" manner, thereby obtaining an error correction code table for the gene sequence. In this way, in the subsequent encoding and decoding process based on the gene sequence code table, encoding and decoding with error correction capabilities can be implemented for the errors that occur.

[0311] Specifically, Figure 8 The figure shows an application example of the method for generating an error correction code table for a gene sequence according to an embodiment of the present application. Figure 8 As shown, when a code table is established by generating a sequence, assuming that n sequences are generated and n=4, a code table that corresponds base segments to 01 segments can be obtained by establishing a correspondence with the binary sequence, where the codeword length L=4, and

[0312] Fig. 9 The figure shows an application example of the coding method based on the gene sequence code table according to the embodiment of the present application. Fig. 9 As shown in the figure, for an input N-bit sequence, by adding a data header and padding zeros at the end, the final sequence length can be divided by n′, where the data header length It represents the rounded value of the logarithm of n' with base 2 as the integer, and the number of zeros padded s satisfies:

[0313] s=(n′-(m+N)mod n′)mod n′

[0314] In this way, the obtained binary sequence is Figure 8 The gene sequence code table shown is encoded bit by bit to obtain an output gene sequence, with a total of N′=(N+s+m) / n′ groups.

[0315] Fig.10 The figure shows an application example of the decoding method based on the gene sequence code table according to the embodiment of the present application. Fig.10 As shown, for the gene sequence to be decoded, it is grouped according to the size of the codeword length L, and the number of groups is equal to the number of encoding groups N′, and if the number of groups is insufficient, the sequence corresponding to 00 (AGCT) is supplemented.

[0316] By comparing block by block, the correct base sequence can be obtained and corresponded to the binary sequence. Then the binary sequences block by block are spliced ​​in order, and the supplementary 0s are removed according to the data header, and the data header is removed at the same time to complete the decoding output.

[0317] Fig.11The figure shows an application example of block-by-block comparison in the decoding method based on the gene sequence code table according to an embodiment of the present application. Fig.11 As shown, the data is input from the lower left and is a set of segments consisting of AGCG (As mentioned above, those skilled in the art can understand that it can also be composed of characters from a given finite character table Σ), let the length of each segment be L (this L is equal to the length of each sequence in the code table), and there are a total of N′ segments (the size of N′ is pre-set before encoding, and the binary sequence is split into N′ parts for encoding. When decoding, it is split according to N′ segments).

[0318] The code table is C.

[0319] The following only discusses ecg's seq1, seq2, mat, startIndex, and endIndex.

[0320] As mentioned above, the input and output of the ecg algorithm are:

[0321] Input the two sequences to be compared seq1 and seq2 respectively (note that the algorithm calculates the edit count that occurs when seq1 becomes seq2, so the order needs to be paid attention to); the starting and ending positions of each comparison startIndex and endIndex; the window width widthWindow to be considered; the depth depthBackTrack to be backtracked; startRow indicates the row from which the dynamic programming matrix will start to be updated; the input dynamic programming matrix mat; chSeq, mask, and shiftAmount obtained from the error correction feature sequence are omitted in the following discussion (this is because this information has been determined during encoding and decoding).

[0322] The output is the updated dynamic programming matrix mat and the initial row number startRow of the next step (startRow will also be omitted in the following text because the current step is always the input of the previous step startRow, and the first step is row 0).

[0323] Therefore, we can only focus on seq1, seq2, mat, startIndex, and endIndex.

[0324]

[0325]

[0326] In the block-by-block alignment algorithm, we first consider the first step. That is, we initialize the dynamic programming matrix mat and make |C| copies, each of which corresponds to a sequence in C. Seq1 takes all the sequences in C, seq2 is s1, startIndex = 0, endIndex = L, as shown in part (1).

[0327] After each codeword C[i] and s1 in the code table are subjected to the ecg algorithm based on the corresponding dynamic programming matrix, the corresponding dynamic programming matrix will be obtained. A total of |C| such matrices will be obtained. We calculate The Loss of , and select the matrix with the largest Loss as the next step, that is:

[0328]

[0329] As shown in part (2).

[0330] The corrected sequence in the first step is C[j1]. (This is equivalent to k=1, because there was no sequence before, so the corrected sequence of group 0 and the sequence of group 0 are both empty)

[0331] Now enter the loop:

[0332] As shown in part (1). Consider the situation at step k. We have the sequence s in the previous step k-1 , the correction sequence C[j k-1 ] and the dynamic programming matrix m from the previous step k-1 . Get C[j k-1 ] and s k-1 Count widthWindow characters from the back to the front as the end of each.

[0333] Next, as shown in part (2), C[j k-1 ] is placed at the end of each sequence of C as the input of seg1, s k-1 Put the end of s k The first part of m is used as the input of seg2. k-1 Make |C| copies, each of which is used as the input to the dynamic programming matrix of the corresponding seg1 and seg2. startIndex is the size of widthWindow, and endIndex is the length of input seg2. In the end, we will get |C| dynamic programming matrices. And select the one with the largest Loss And output the corresponding dynamic programming matrix This part is shown in the gray block in the middle of the figure. (Note: The purpose of adding the end of the previous step in this step is to count the errors that occur at the junction of two sequences. When the length of one of the feature sequences in the error correction feature sequence exceeds 1, it is possible that the error occurs at the junction, that is, one part is at the end of the kth sequence and the rest is at the beginning of the k+1th sequence. The range of "spanning" errors that we need to pay attention to is exactly the size of widthWindow.).

[0334] Then, as shown in part (3), as shown in the gray block in the upper right corner, we get the sequence s at step k k , the correction sequence C[j k ], the dynamic programming matrix m of step k k .

[0335] Finally, as shown in part (4), we need to first determine Loss(m k ) is 0:

[0336] If it is 0, the base sequence corresponding to the all-0 binary sequence is supplemented, and then all previously corrected fragments are output.

[0337] If it is not zero, it is necessary to determine whether the correction has been completed (i.e., is the currently completed sequence the last sequence, k == N′?):

[0338] If they are equal, the task is completed and the result is output;

[0339] Otherwise, return to the upper left corner and proceed to the next step of calculation.

[0340] Two error correction examples of the encoding and decoding method based on the gene sequence code table according to the embodiment of the present application are described below. In the first example, the error correction feature sequence is as follows:

[0341] ((((″,′*′),),1),(((′*′,″),),1))

[0342] This means that the obtained code table can correct 1 insertion-1 deletion. Let the size of the generated code table be 8 and the run limit be 2, and obtain the following results:

[0343] ATCGACACAGCT GTGCTGACAGAT CGTGGATCGAAT TGTCCATCCGAT AAGCTCACTGCT GACCGTTCAGAT CCACCAGACAAT TGGAGTTCAGCT

[0344] The order is from left to right and from top to bottom. If you want to encode 001 101 111, then there is a sequence

[0345] GTGCTGACAGAT GACCGTTCAGAT TGGAGTTCAGCT

[0346] If the following error occurs:

[0347] GTG_TGACAGAT GACCGTTCAGAT TG C GAGTTCAGCT

[0348] According to the decoding algorithm, the loss function of each group is as follows:

[0349] [0.0,0.6,0.0,0.0,0.0,0.0,0.0,0.0]

[0350] [0.0,0.0,0.0,0.0,0.0,0.6,0.0,0.0]

[0351] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.2]

[0352] According to the non-zero indexes, they are 1 5 7, and the corresponding binary numbers are 001 101 111, which is the original input.

[0353] If more errors occur, such as:

[0354] GTG_TGA_AGAT GACCGTTCAGAT TG C GAGTTCAGCT

[0355] The loss function for each step of the decoding algorithm is as follows:

[0356] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]

[0357] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]

[0358] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]

[0359] Therefore, if the number of errors exceeds the preset value, binary 000 000 000 is returned.

[0360] In the second example, let the error correction feature sequence be as follows:

[0361] ((((′A′,″),),1),(((′G′,′C′),),1))

[0362] This means that the encoding completed by this code table can completely correct one deletion of 'A' and one replacement of 'G' with 'C'. The code table size and run length restrictions are the same as above, and the following results are obtained:

[0363] ATGCGGTA GAGGCTTA CACCGTTA TTCGAGCA ACGAGTTC GACGTCTA CTCCTGTA TACCAGGA

[0364] Same order as above.

[0365] Still encoding the binary 001 101 111, there is a sequence:

[0366] GAGGCTTA GACGTCTA TACCAGGA

[0367] If there is one deletion of 'A' and one replacement of 'G' with 'C', as follows:

[0368] GAG C CTTA G_CGTCTA TACCAGGA

[0369] The loss function for each group is:

[0370] [0.0,0.6,0.0,0.0,0.0,0.0,0.0,0.0]

[0371] [0.0,0.0,0.0,0.0,0.0,0.2,0.0,0.0]

[0372] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.2]

[0373] The corresponding binary is 001 101 111.

[0374] If another error occurs, such as inserting a 'C'

[0375] GAGGCTTAGA C CGTCTATACCAGGA

[0376] The loss function for each group is:

[0377] [0.0,1.4,0.0,0.0,0.0,0.0,0.0,0.0]

[0378] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]

[0379] [0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0]

[0380] Therefore, the value returned is 001 000 000, and the number of errors exceeds the preset value.

[0381] In summary, according to the embodiment of the present application, the error correction code table generation method for gene sequences and the encoding and decoding method based on the gene sequence code table can be based on specific error correction requirements, that is, the error correction error type and the number of errors corresponding to each type, to achieve code table generation with an appropriate size, and encoding and decoding based on the generated code table, thereby achieving high-performance code table generation and encoding and decoding processing.

[0382] The basic principles of the present application are described above in conjunction with specific embodiments. However, it should be noted that the advantages, strengths, effects, etc. mentioned in the present application are only examples and not limitations, and it cannot be considered that these advantages, strengths, effects, etc. are required by each embodiment of the present application. In addition, the specific details disclosed above are only for the purpose of illustration and ease of understanding, not for limitation, and the above details do not limit the present application to being implemented by adopting the above specific details.

[0383] The block diagrams of the devices, apparatuses, equipment, and systems involved in this application are only illustrative examples and are not intended to require or imply that they must be connected, arranged, and configured in the manner shown in the block diagram. As will be appreciated by those skilled in the art, these devices, apparatuses, equipment, and systems can be connected, arranged, and configured in any manner. Words such as "including", "comprising", "having", etc. are open words, referring to "including but not limited to", and can be used interchangeably with them. The words "or" and "and" used here refer to the words "and / or" and can be used interchangeably with them, unless the context clearly indicates otherwise. The words "such as" used here refer to the phrase "such as but not limited to", and can be used interchangeably with them.

[0384] It should also be noted that in the apparatus, device and method of the present application, each component or each step can be decomposed and / or recombined. Such decomposition and / or recombination should be regarded as equivalent solutions of the present application.

[0385] The above description of the disclosed aspects is provided to enable any person skilled in the art to make or use the present application. Various modifications to these aspects will be readily apparent to those skilled in the art, and the general principles defined herein may be applied to other aspects without departing from the scope of the present application. Therefore, the present application is not intended to be limited to the aspects shown herein, but rather to the widest scope consistent with the principles and novel features disclosed herein.

[0386] The above description has been given for the purpose of illustration and description. In addition, this description is not intended to limit the embodiments of the present application to the forms disclosed herein. Although multiple example aspects and embodiments have been discussed above, those skilled in the art will recognize certain variations, modifications, changes, additions and sub-combinations thereof.

Claims

1. A method for generating an error correction code table, comprising: Determine a finite alphabet of a sequence code as an error correction target and a code table size of a corresponding error correction code, and generate an error correction feature sequence based on the finite alphabet and the error correction error type, wherein the error correction feature sequence is used to define the number of error correction errors of the code table based on the error correction error type; Determine a feature sequence, a feature mask and a shift amount based on the error correction feature sequence; Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount; and In response to the loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, two code words in a code table corresponding to each dynamic programming matrix completed by iteration are determined to obtain an error correction code table.

2. The method for generating an error correction code table according to claim 1, wherein: Generating an error correction feature sequence based on the finite alphabet and the error correction error type includes: generating a complete alphabet based on the limited alphabet; Defining a feature string space by a disjoint union of multiple complete character tables, wherein the feature string space includes an empty string representation; The characteristic sequence pair representing the error type of error correction is defined by two different elements in the monoid defined by character concatenation on the characteristic string space; Defining an error correction feature unit based on a plurality of said feature sequence pairs and an allowable upper limit as a positive integer greater than or equal to zero; and An error correction feature sequence is defined based on a sequence consisting of multiple error correction feature units.

3. The method for generating an error correction code table according to claim 2, wherein: The characteristic sequence pair indicating the type of error correction error includes: represents the replacement (p1, p2), p i ≠-,i=1,2; represents the insertion of (-,p2), p2≠-; and Represents deleted (p1,-), p1≠-; Where p1, p2 are two different elements in the monoid, and - is the empty string.

4. The method for generating an error correction code table according to claim 1, wherein: After the error correction feature sequence is generated, the method further includes performing symmetric processing on the error correction feature sequence, including: For the error correction feature sequence, start traversing each error correction feature unit; For each error correction feature unit, multiply the allowable upper limit by two; and For each feature sequence pair in the error correction feature unit, determine whether a symmetric feature sequence pair that is arranged opposite to two elements of the feature sequence is in the error correction feature unit, and in response to the symmetric feature sequence pair not being in the error correction feature unit, place the symmetric feature sequence pair in the error correction feature unit.

5. The method for generating an error correction code table according to claim 1, wherein: Determining a feature sequence based on the error correction feature sequence includes: The feature sequence pair corresponding to the unit feature in each error correction feature unit in the error correction feature sequence is projected into the feature sequence.

6. The method for generating an error correction code table according to claim 1, wherein: Determining the feature mask and the translation amount based on the error correction feature sequence includes: An array for determining an allowable upper limit based on the error correction feature sequence; Determine an edit information tuple based on the array of allowable upper limits; Determining a bit array index mapping based on the arrays of the edit information tuple and the allowable upper limit, the bit array index mapping being used to map the edit information tuple into a zero-based order based on a lexicographic order; and The feature mask and the translation amount are determined based on the bit array index mapping.

7. The method for generating an error correction code table according to claim 1, wherein: Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation comprises: When the size of the code table is n×l, n empty character strings are initialized as code words of the code table, and the length of each code word is l; Initialize n(n-1) / 2 dynamic programming matrices, wherein each dynamic programming matrix corresponds to two codewords of the code table, and its size is (backtracking depth+1, 2*window size+1), each position of the dynamic programming matrix is ​​a bit array, and the 0th position of the bit array of the 0th row and the middle column of the initialized dynamic programming matrix is ​​1; and Initialize n(n-1) / 2 starting rows.

8. The method for generating an error correction code table according to claim 7, wherein: Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation comprises: The dynamic programming matrix is ​​iterated using the edit computation graph algorithm, where in the i-th iteration: For the starting row number determined by the starting coordinates and the ending coordinates of the codeword comparison, determine the comparison index of the two codewords based on i and the backtracking depth, and determine a comparison string pair of the two codewords based on the comparison index; In response to two compared string pairs being equal, updating the dynamic programming matrix using a first matrix updating algorithm; and In response to the two compared character string pairs being unequal and conforming to the feature sequence pair, the dynamic programming matrix is ​​updated using a second matrix updating algorithm according to the feature mask and the translation amount.

9. The method for generating an error correction code table according to claim 8, wherein: Updating the dynamic programming matrix using a first matrix update algorithm includes: Determine the coordinates of the dynamic programming matrix by subtracting one from the comparison indexes of the two codewords; Extract the data of the dynamic programming matrix at the coordinates by masking; and Shift to the higher bit and put it into the position of the dynamic programming matrix with the comparison index of the two code words as the coordinate through OR operation.

10. The method for generating an error correction code table according to claim 8, wherein: Updating the dynamic programming matrix using a second matrix update algorithm according to the feature mask and the translation includes: Determine the coordinates of the dynamic programming matrix by subtracting the current feature sequence pair from the comparison index of the two code words; Extract the data of the dynamic programming matrix at the coordinates by using a mask; an amount by which the translation amount is translated toward the upper position; and The codewords are placed in the dynamic programming matrix at a position with the comparison index of the two codewords as coordinates through an OR operation.

11. The method for generating an error correction code table according to claim 10, wherein: The loss function of the dynamic programming matrix is ​​expressed as: Wherein d is the backtracking depth, w is the window size, s is the starting function, and (i, j) is the comparison index of the two codewords.

12. The method for generating an error correction code table according to claim 1, further comprising: When the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, the starting coordinates and the ending coordinates of the dynamic programming matrix iteration are determined based on the length of each codeword in the new code table and the window width; as well as A new starting row number is determined based on the starting coordinate and the ending coordinate, and a dynamic programming matrix is ​​continuously iterated based on the dynamic programming matrix corresponding to the generated code table based on the new starting row number to obtain the new code table.

13. A coding method based on a generated code table, comprising: Determine a floor value of a logarithmic function with base 2 of a code table size of the generated code table; Preprocess the input binary sequence to pad the end with zeros and record the number of padding zeros; and The pre-processed input binary sequence is grouped and encoded with the number of the rounded-down values ​​as a group.

14. A decoding method based on a generated code table, comprising: Dividing the received sequence to be decoded into a predetermined number of groups, where the predetermined number is the number of coding groups, and the length of each group is the length of each codeword in the generated code table; Starting from the first group, each group is compared with the code table to determine the code word in the code table corresponding to the group; For the kth group in the sequence to be decoded and the corresponding error correction sequence index is i, the codeword in the code table is added with the prefix of the codeword minus the window size, and matched with the sequence of the previous correction sequence minus the start of the window size and the length of the codeword multiplied by k; and When all groups of the sequence to be decoded are matched, all correct sequences are concatenated to obtain a decoded sequence.

15. The decoding method based on the generated code table as claimed in claim 14, further comprising over-range processing: When the character string corresponding to the group is all 0, the last zero sequence is taken to correspond to it, and it is subsequently defaulted to the first code word in the code table; When there are multiple non-zero characters in the character string corresponding to the group, the sequence with the largest loss function is selected to correspond to it.

16. A method for generating an error correction code table for a gene sequence, comprising: Determine the size of the error correction code table for the gene sequence, and generate an error correction feature sequence based on the gene sequence {A, T, C, G} as a finite character table and the error correction error types and error correction error numbers corresponding to the error correction target of the gene sequence; Determine a feature sequence, a feature mask and a shift amount based on the error correction feature sequence; Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation amount; and In response to the loss function of the dynamic programming matrix being less than or equal to a predetermined threshold, two code words in a code table corresponding to each dynamic programming matrix completed by iteration are determined to obtain an error correction code table for a gene sequence.

17. The method for generating an error correction code table for a gene sequence according to claim 16, wherein: Generating an error correction feature sequence based on the gene sequence {A, T, C, G} as a finite character table and the error correction error types and error correction error numbers corresponding to the gene sequence error correction target includes: Generate a complete alphabet based on the gene sequence {A, T, C, G} as a finite alphabet; Defining a feature string space by a disjoint union of multiple complete character tables, wherein the feature string space includes an empty string representation; The characteristic sequence pair representing the error type of error correction is defined by two different elements in the monoid defined by character concatenation on the characteristic string space; Defining an error correction feature unit based on a plurality of said feature sequence pairs and an allowable upper limit as a positive integer greater than or equal to zero; and An error correction feature sequence is defined based on a sequence consisting of multiple error correction feature units.

18. The method for generating an error correction code table for a gene sequence according to claim 17, wherein: The characteristic sequence pair indicating the type of error correction error includes: represents the replacement (p1, p2), p i ≠-,i=1,2; represents the insertion of (-,p2), p2≠-; and Represents deleted (p1,-), p1≠-; Where p1, p2 are two different elements in the monoid, and - is the empty string.

19. The method for generating an error correction code table for a gene sequence according to claim 17, wherein: After the error correction feature sequence is generated, the method further includes performing symmetric processing on the error correction feature sequence, including: For the error correction feature sequence, start traversing each error correction feature unit; For each error correction feature unit, multiply its allowable upper limit by two; and For each feature sequence pair in the error correction feature unit, determine whether a symmetric feature sequence pair that is arranged opposite to two elements of the feature sequence is in the error correction feature unit, and in response to the symmetric feature sequence pair not being in the error correction feature unit, place the symmetric feature sequence pair in the error correction feature unit.

20. The method for generating an error correction code table for a gene sequence according to claim 16, wherein: Determining a feature sequence based on the error correction feature sequence includes: The feature sequence pair corresponding to the unit feature in each error correction feature unit in the error correction feature sequence is projected into the feature sequence.

21. The method for generating an error correction code table for a gene sequence according to claim 16, wherein: Determining the feature mask and the translation amount based on the error correction feature sequence includes: An array for determining an allowable upper limit based on the error correction feature sequence; Determine an edit information tuple based on the array of allowable upper limits; determining a bit array index mapping based on the arrays of the edit information tuple and the permissible upper limit, the bit array index mapping being used to map the edit information tuple into a zero-based order based on a lexicographic order; and, The feature mask and the translation amount are determined based on the bit array index mapping.

22. The method for generating an error correction code table for a gene sequence according to claim 16, wherein: Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation comprises: When the size of the code table is n×l, n empty character strings are initialized as code words of the code table, and the length of each code word is l; Initialize n(n-1) / 2 dynamic programming matrices, wherein each dynamic programming matrix corresponds to two codewords of the code table, and its size is (backtracking depth+1, 2*window size+1), each position of the dynamic programming matrix is ​​a bit array, and the 0th position of the bit array of the 0th row and the middle column of the initialized dynamic programming matrix is ​​1; and, Initialize n(n-1) / 2 starting rows.

23. The method for generating an error correction code table for a gene sequence according to claim 21, wherein: Iterating a dynamic programming matrix based on the feature sequence, the feature mask and the translation comprises: The dynamic programming matrix is ​​iterated using the edit computation graph algorithm, where in the i-th iteration: For the starting row number determined by the starting coordinates and the ending coordinates of the codeword comparison, determine the comparison index of the two codewords based on i and the backtracking depth, and determine a comparison string pair of the two codewords based on the comparison index; In response to two compared string pairs being equal, updating the dynamic programming matrix using a first matrix updating algorithm; and, In response to the two compared character string pairs being unequal and conforming to the feature sequence pair, the dynamic programming matrix is ​​updated using a second matrix updating algorithm according to the feature mask and the translation amount.

24. The method for generating an error correction code table for a gene sequence according to claim 23, wherein: Updating the dynamic programming matrix using a first matrix update algorithm includes: Determine the coordinates of the dynamic programming matrix by subtracting one from the comparison indexes of the two codewords; Extract the data of the dynamic programming matrix at the coordinates by masking; and, Shift to the higher bit and put it into the position of the dynamic programming matrix with the comparison index of the two code words as the coordinate through OR operation.

25. The method for generating an error correction code table for a gene sequence according to claim 22, wherein: Updating the dynamic programming matrix using a second matrix update algorithm according to the feature mask and the translation includes: Determine the coordinates of the dynamic programming matrix by subtracting the current feature sequence pair from the comparison index of the two code words; Extract the data of the dynamic programming matrix at the coordinates by using a mask; an amount by which the translation amount is translated toward the upper position; and, The codewords are placed in the dynamic programming matrix at a position with the comparison index of the two codewords as coordinates through an OR operation.

26. The method for generating an error correction code table for a gene sequence according to claim 16, wherein: The loss function of the dynamic programming matrix is ​​expressed as: Wherein d is the backtracking depth, w is the window size, s is the starting function, and (i, j) is the comparison index of the two codewords.

27. The method for generating an error correction code table for a gene sequence according to claim 16, further comprising: When the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, the starting coordinates and the ending coordinates of the dynamic programming matrix iteration are determined based on the length of each codeword in the new code table and the window width; and A new starting row number is determined based on the starting coordinate and the ending coordinate, and a dynamic programming matrix is ​​continuously iterated based on the dynamic programming matrix corresponding to the generated code table based on the new starting row number to obtain the new code table.

28. The method for generating an error correction code table for a gene sequence according to claim 27, wherein: Before determining the starting coordinates and the ending coordinates of the dynamic programming matrix iteration based on the length of each codeword and the window width in the new code table, the method further includes: When the length of each codeword in the generated code table is L and a new code table with a codeword length greater than L is to be generated, a base sequence is added as a suffix to each generated codeword according to a predetermined probability corresponding to each gene base based on an incremental strategy.

29. The method for generating an error correction code table for a gene sequence as claimed in claim 28, wherein: The increment strategy is determined based on a suffix increment constraint corresponding to the gene sequence, wherein the suffix increment constraint includes a GC ratio constraint and a continuous identical base constraint.

30. The method for generating an error correction code table for a gene sequence according to claim 28, wherein: The incremental strategy is based on a Monte Carlo strategy, in which the probability of each base is 0.25 when adding a suffix in each round, and the probability is incrementally constrained based on a loss function according to the codewords of the generated code table.

31. The method for generating an error correction code table for a gene sequence as claimed in claim 28, wherein: The incremental strategy is based on a deep reinforcement learning strategy, and the deep reinforcement learning strategy includes: Inputting the code table into a deep neural network to obtain a feature matrix, wherein the feature matrix represents the probability of adding different bases to each code word of the code table; Obtaining a result matrix from the feature matrix based on a suffix increment constraint; Generate candidate suffixes based on the result matrix, and calculate the optimal suffix, reward and action by minimizing the loss function; and The deep reinforcement learning strategy is optimized through policy gradient according to the reward and action.

32. The method for generating an error correction code table for a gene sequence as claimed in claim 31, wherein: The rewards are: Obtain the code table, code table length L0, and loss function R0 through the Monte Carlo method; In the reinforcement learning process, the new loss function l after adding the suffix in the i-th step is compared with the loss function of the corresponding step of R0 to determine the reward: In response to l≤L0, the reward is (R0[i]-l) / (R0[i]); In response to l>L0, the reward is -1.

33. A coding method based on a gene sequence code table, comprising: Determine the rounded-down value of a logarithmic function with base 2 of the code table size of the gene sequence code table; Preprocess the input binary sequence to pad the end with zeros and record the number of padding zeros; and The pre-processed input binary sequence is grouped and encoded with the number of the rounded-down values ​​as a group.

34. The encoding method based on the gene sequence code table as claimed in claim 33, wherein: When the size of the code table of the generated code table is n×l and the length of the binary sequence to be encoded is N, N satisfies: in And s=(n′-(m+N)mod n′)mod n′.

35. A decoding method based on a gene sequence code table, comprising: Dividing the received gene sequence to be decoded into a predetermined number of groups, the predetermined number being the number of coding groups, and the length of each group being the length of each codeword in the gene sequence code table; Starting from the first group, each group is compared with the gene sequence code table to determine the code word corresponding to the group in the gene sequence code table; For the kth group of the gene sequence to be decoded and the corresponding error correction sequence index is i, the codeword in the gene sequence code table is added with the prefix of the codeword minus the window size, and matched with the sequence of the previous correction sequence minus the start of the window size and the length of the codeword multiplied by k; and When all groups of the gene sequences to be decoded are matched, all correct sequences are spliced ​​to obtain a decoded gene sequence.

36. The decoding method based on the generated code table as claimed in claim 35, further comprising over-range processing: When all the character strings corresponding to the group are 0, the last zero sequence is taken to correspond to it, and is subsequently defaulted to be the first codeword in the code table; and When there are multiple non-zero characters in the character string corresponding to the group, the sequence with the largest loss function is selected to correspond to it.

Images