Data processing method for determining DNA (Deoxyribose Nucleic Acid) organism through biological nanopores

Abstract

The invention discloses a data processing method for determining a DNA nucleic acid organism through a biological nanopore, which comprises a series of key steps of data correction, noise filtering, event screening, level searching, sequence alignment and the like from the acquisition of original current data to the generation of a final DNA base sequence and the determination of the organism, and constructs a set of complete data processing flow. According to the method, comprehensive processing and analysis of biological nanopore current signals are realized, and DNA sequence information can be acquired in real time by adopting a mode of simultaneously acquiring and analyzing, so that the detection efficiency is greatly improved. Analysis can be carried out without waiting for the completion of all data acquisition, so that a lot of time is saved; the biological nanopore is used as a detection means, so that the sample preparation requirement is low, complicated sample treatment is not needed, and the data analysis process is clear.

Description

Technical Field

[0001] The present invention belongs to the technical field of biological information nanopore sequencing, and in particular relates to a data processing method for determining a DNA nucleic acid organism using a biological nanopore. Background Art

[0002] Bio-nanopore is a self-assembled protein pore. This protein pore is embedded in a phospholipid bilayer membrane or a polymer membrane. The current can be detected by adding an ion buffer and applying a suitable voltage. When the substance to be detected is added to the system, it can be detected that a blocking current is generated when the substance to be detected passes through. In theory, the blocking current data is the current signal generated by the substance to be detected. Other abnormal events caused by the environment or experiment cannot be ruled out, but generally this situation rarely occurs.

[0003] At present, the principle of conducting experiments and tests is well known in the industry, but only the UK and the company Oxford Nanopore are proficient in and commercially apply this technology. At present, in China, this principle still has a lot of room for exploration in experimental operation technology, application field development and experimental data processing. The maturity of experimental operation technology, application field development and experimental data processing technology will promote the development of related applications of gene sequencing in China and even internationally.

[0004] Bio-nanopore technology can be used to sequence DNA quickly and accurately, so it can accelerate the DNA sequencing process in genomic research and help scientists better understand genome structure, function and variation. In medical diagnosis, this method can be used to detect gene mutations, diagnose genetic diseases and personalized treatment. Bio-nanopore technology can help scientists better understand the interaction between drugs and target DNA, thereby accelerating the drug development process. In addition, the technology can also be used to detect the purity and quality of drugs and monitor DNA impurities in the biopharmaceutical process. Bio-nanopore technology can be used to detect microbial DNA in the environment, thereby helping to monitor environmental pollution, predict disease outbreaks and assess the health of ecosystems. In terms of food safety, the technology can be used to detect pathogenic microorganisms or genetically modified ingredients in food to ensure the quality and safety of food. With the advancement of bio-nanopore technology, the cost of personal genome sequencing is gradually decreasing, which will provide more possibilities for personalized health management. Individuals can understand their genetic characteristics, disease risks and drug reactions by regularly sequencing their own genomes, so as to take corresponding preventive and treatment measures to improve their health level and quality of life. Bio-nanopore technology can be used in the field of forensics to help detect and identify DNA evidence left at crime scenes, thereby helping to solve criminal cases and ensure judicial justice. Compared with traditional DNA analysis technology, biological nanopore technology has higher resolution and sensitivity and can provide more reliable evidence. With the development of application technology, the application scope of biological nanopore for DNA nucleic acid detection of organisms is becoming more and more extensive, but the data processing algorithm or implementation plan that comes with the technology is still unclear, and data processing methods are urgently needed to accelerate the advancement of biological nanopore detection technology for organisms. Summary of the invention

[0005] In view of the defects and problems existing in the existing DNA nucleic acid organism detection, the present invention provides a data processing method for determining DNA nucleic acid organisms using biological nanopores, aiming to improve sequencing accuracy and reliability, simplify the data analysis process, and accelerate the optimization and promotion of nanopore sequencing technology.

[0006] The solution adopted by the present invention to solve the technical problem is: a data processing method for determining DNA nucleic acid organisms using a biological nanopore, comprising the following steps:

[0007] Step 1: According to the characteristics of the biological nanopore, the DNA sequence to be detected is divided into a combination of several bases, and then the DNA sequence to be detected is subjected to biological nanopore detection; and the original current change data of the DNA sequence to be detected is obtained;

[0008] Step 2, correcting the invalid values ​​in the original current change data by using a correction algorithm;

[0009] Step 3, performing high-frequency noise filtering on the corrected current change data;

[0010] Step 4: divide the processed current change data into opening current and blocking current, and filter the events that meet the conditions;

[0011] Step 5: perform level search and fitting for the events that meet the conditions. First, identify the level, determine the boundary of the current resistance level generated by different base combinations, and then remove the false level.

[0012] Step 6: Establish a reference level database, select DNA samples with known sequences as test sets, perform nanopore sequencing and level extraction on them, compare the level data of the test set samples with the reference level database, obtain the annotated DNA sequence, compare the annotated DNA sequence with the known actual sequence of the test set samples, calculate the sequence similarity, analyze the type and source of the mismatched bases, evaluate the accuracy of the reference database, and optimize and improve the reference database according to the evaluation results, such as removing erroneous data and adding new samples;

[0013] Step seven, based on the Needleman-Wunsch and Smith-Waterman sequence alignment algorithms, compare the level data obtained from the DNA sequencing to be tested in step five with the level data in the reference level database, find the most similar sequence, infer the DNA sequence of the DNA to be tested, add the level of the newly measured DNA sequence to the reference level database, and recalculate the median of different sequence combinations to update and improve the reference level database and improve its accuracy.

[0014] In the step 1, the DNA sequence to be detected is detected by biological nanopore using a sampling frequency of 500k.

[0015] In the step 2, the invalid value is a value that exceeds a preset percentage of the normal value, and the correction algorithm includes one of a median filter and a weighted average algorithm.

[0016] In step 3, a Bayesian low-pass filter is used to filter the high-frequency noise of the corrected original current change data. The sampling frequency of this project is 500KHz, the filter adopts Bessel 4th order, and the cut-off frequency is 50KHz (ie Fs / 10).

[0017] In step 4, event refers to a blocking current event generated when DNA passes through the pore, and the screening of event includes the following steps:

[0018] Step 1: First, based on a broad threshold detection, select the extracted event with a duration Time>1S and a current range between (0.1, 0.75) of the detected opening current;

[0019] Step 2, initializing the local opening current value (L initialization local) according to the number of events detected (L) and the opening current detected;

[0020] Step 3, then according to each local open state (L i ) and each Event(L i ) to analyze each event one by one and further screen suitable events.

[0021] The event further screening includes the following steps:

[0022] Step 1: Check the states before and after each numbered event in event: event_before and event_after, which represent the open states before and after the event respectively, and the initial value is NaN.

[0023] Step 2: If event is not the first state, the median of the open state before event is calculated as event_before; if event is not the last state, the median of the state after event is calculated as event_after.

[0024] Step 3: If the current event number is greater than 2, and there is data in the first two open states, try to use the exponential decay function to fit the data of the previous state and update it to event_before.

[0025] Step 4: If the current event number is less than the open state number and there is data in the current state, try to use the exponential decay function to fit the data of the next state and update it to event_after.

[0026] Step 5, exception handling, if event_before and event_after abnormally deviate from the 3sigma value of the open-state current detected by event, they are also set to NaN.

[0027] The establishment of the reference level database in step 6 includes the following steps:

[0028] Step 1, select a DNA sample with a known sequence for nanopore sequencing to obtain raw current data;

[0029] Step 2: Process and analyze the current data to extract reliable level information;

[0030] Step 3, repeat 1-2 multiple times, accumulate level data, and establish a reference level database.

[0031] Step 4: Associate the levels in the reference level database with their corresponding known DNA sequences to form a "level-sequence" mapping relationship.

[0032] Beneficial effects of the present invention: The data processing method for determining DNA nucleic acid organisms by a biological nanopore provided by the present invention has the following beneficial effects:

[0033] 1. From the collection of raw current data to the generation of the final DNA base sequence, it covers a series of key steps such as data correction, noise filtering, event screening, level search, sequence alignment, etc., and builds a complete data processing flow to achieve comprehensive processing and analysis of biological nanopore current signals.

[0034] 2. The method of collecting and analyzing at the same time can obtain DNA sequence information in real time, which greatly improves the detection efficiency. There is no need to wait until all data collection is completed before analysis, which saves a lot of time.

[0035] 3. Using biological nanopores as the detection method, the sample preparation requirements are low, no complicated sample processing is required, and the data analysis process is clear and straightforward.

[0036] 4. By setting the threshold reasonably and using statistical methods to identify valid events, the interference of false events is eliminated. Level fitting fully considers the frameshift situation, allowing jumps and retreats during sequence alignment, and improving the accuracy of sequence generation.

[0037] 5. The database is continuously improved through level extraction of known sequence samples, and the newly measured levels are continuously supplemented and updated, forming a dynamic association between sequence and level, providing a reliable reference for sequence inference of unknown DNA samples. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] Figure 1 It is a schematic diagram of current change data preprocessing of the present invention.

[0039] Figure 2 It is a schematic diagram of the classification of the opening current and the blocking current of the present invention.

[0040] Figure 3 It is a schematic diagram of event generation of the present invention.

[0041] Figure 4 It is a schematic diagram of data window traversal of the present invention.

[0042] Figure 5It is a schematic diagram of the recursive search for the optimal t2 segmentation point of the present invention.

[0043] Figure 6 It is a schematic diagram of whether the window of the present invention has transition.

[0044] Figure 7 It is a schematic diagram of the optimal transition position of the present invention.

[0045] Figure 8 It is a schematic diagram of the Level fitting result of the present invention.

[0046] Fig. 9 It is the best comparison logic schematic diagram of the present invention.

[0047] Fig.10 It is the optimal route schematic diagram of the present invention. DETAILED DESCRIPTION

[0048] The present invention is further described below in conjunction with the accompanying drawings and embodiments.

[0049] Example 1

[0050] This embodiment provides a data processing method for determining a DNA nucleic acid organism using a biological nanopore, and other features include the following steps:

[0051] Step 1: According to the characteristics of the biological nanopore, the DNA sequence to be detected is divided into a combination of several bases. For example, the DNA sequence to be detected can be designed as a group of 4 bases, a group of 5 bases, or a group of 6 bases, etc., so the sequence combinations are 44, 45, 46, etc.; then the DNA sequence to be detected is subjected to biological nanopore detection; the original current change data of the DNA sequence to be detected is obtained; the DNA sequence to be detected is subjected to biological nanopore detection using a sampling frequency of 500k.

[0052] Step 2: Correct the invalid values ​​in the original current change data through a correction algorithm; an invalid value is a value that exceeds a preset percentage of a normal value, for example, a value that exceeds a preset 30% of a normal value is an invalid value, and the correction algorithm includes one of a median filter and a weighted average algorithm.

[0053] Step 3: Perform high-frequency noise filtering on the corrected current change data, such as Figure 1 As shown, the Bayesian low-pass filter is used to filter the high-frequency noise of the corrected raw current change data. The sampling frequency of this project is 500KHz, the filter adopts Bessel 4th order, and the cut-off frequency is 50KHz (i.e. Fs / 10).

[0054] Step 4: Figure 2As shown in the figure, the processed current change data is divided into the opening current and the blocking current. The states at this time are called the open state and the blocking state, and the events that meet the conditions are screened; the event refers to the blocking current event generated when there is DNA passing through the pore, such as Figure 3 The figure shows the generation of an event. The event screening includes the following steps:

[0055] Step 1: First, based on a broad threshold detection, select the extracted event with a duration Time>1S and a current range between (0.1, 0.75) of the detected opening current;

[0056] Step 2, initializing the local opening current value (L initialization local) according to the number of events detected (L) and the opening current detected;

[0057] Step 3, then according to each local open state (L i ) and each Event(L i ) to analyze each event one by one and further screen suitable events.

[0058] Further event screening includes the following steps:

[0059] Step 1: Check the states before and after each numbered event in event: event_before and event_after, which represent the open states before and after the event respectively, and the initial value is NaN.

[0060] Step 2: If event is not the first state, the median of the open state before event is calculated as event_before; if event is not the last state, the median of the state after event is calculated as event_after.

[0061] Step 3: If the current event number is greater than 2, and there is data in the first two open states, try to use the exponential decay function to fit the data of the previous state and update it to event_before.

[0062] Step 4: If the current event number is less than the open state number and there is data in the current state, try to use the exponential decay function to fit the data of the next state and update it to event_after.

[0063] Step 5, exception handling, if event_before and event_after abnormally deviate from the 3sigma value of the open-state current detected by event, they are also set to NaN.

[0064] Step 5: perform level search and fitting for the event that meets the conditions. The level search includes the following steps:

[0065] Step 1: Perform a level search for valid event data; first define the data window Traversal The entire event, such as Figure 4 shown.

[0066] Step 2, define t1 and t3 as window time positions, and find the optimal t2 split point through recursive search, such as Figure 5 As shown, the evaluation condition is that the smaller the better, see formula 3 in the manual for details.

[0067] Step 3, then search for t2 in the window (t1, t2) and window (t2, t3) respectively. When there is no transition in the window, expand the search range; when there are multiple transitions in the window, narrow the search range, such as Figure 6 shown.

[0068] Step 4: Finally, by repeating the second and third steps, all optimal transition positions are found. The result is as follows: Figure 7 , you can search for all levels in the event. Of course, when evaluating the optimal transition position, the conditions are as follows: P(t1, t2, t3) is less than the set threshold, see formula 3 in the manual for details.

[0069] The detailed methods and formulas involved in level search mainly involve two steps. One is to identify the level and determine the boundary within the time trajectory of the level. The second is to remove false levels. The specific steps are as follows:

[0070] Step 1, determine the location of the level, first check a part of the current trajectory, and divide this part of the current into two parts. Assuming that the sampled current in each level follows a Gaussian distribution, calculate the total probability that the two currents come from two different Gaussian distributions, and divide this total probability by the probability of the null hypothesis that the combination of the two parts comes from a single Gaussian. For ease of calculation, use the logarithmic probability. For the current intercepted in a given window, the average value of the observed current is I mean , width is w, for this part (T time) the single current I T , at time t, the logarithmic probability (density) is:

[0071]

[0072] Step 2: Find the observations between T1 and T2 that have I mean The total log probability of a single level defined by and w, summing the formula, that is, Formula 1 between T1 and T2, and using sigma to define const, gives:

[0073] LogP(I t1，t2 )=(t2-t1)Logo(t2,t1)+const Formula 2

[0074] in:

[0075]

[0076] Step 3, calculate the total log probability, one of the two parts between T1 and T2 and between T2 and T3. To compare the log probabilities, subtract the total log probability of the null hypothesis by calculating the equation. Combining all the probabilities produces a comparison measure and represents it as LogP(T1,T2,T3):

[0077] LogP(t1,t2,t3)=(t2-t1)Logσ(t2,t1)+(t3-t2)Logσ(t3,t2)-(t3-t1)Logσ(t3,t1) Formula 3

[0078] in:

[0079]

[0080] The minimum LogP(t1, t2, t3) produced by T2 indicates the location where an energy level transition may occur between the currently observed t1 and t3. In the level search, start from the given time window [T1, T3] and search for T2 that minimizes LogP. If min(LogP) is less than the specified threshold (LogP = -50), there is a level transition at T2. If min(LogP) is higher than the threshold of the original time window, there is no transition between t1 and t3. Consider increasing the value of t3 to add a larger window.

[0081] Step 4, delete false levels. Not all levels meet the requirements. Due to other reasons, such as: false levels caused by DNA blocking, environmental disturbances, etc., delete ① the duration is less than 500uS, ① the average value exceeds the threshold range [9-72pA] or the average value error is greater than 5pA. Finally, use the logarithm calculated by the formula. In addition, similar levels will be identified in formula 3 and the levels of the levels will be merged.

[0082] Step 6: Establish a reference level database. After determining the event of the DNA to be detected, search for the level of the event. The level is generated because of different base combinations. For example, each time the biological nanopore reads 4 bases, the levels generated by different combinations of 256 bases such as AAAA, AAAC, AAAG, etc. are basically different, which is mainly reflected in the amplitude and duration of the current blockade. Figure 8 As shown, a standard level reference database is constructed. DNA samples with known sequences are selected as test sets, and nanopore sequencing and level extraction are performed on them. The level data of the test set samples are compared with the reference level database to obtain annotated DNA sequences. The annotated DNA sequences are compared with the known actual sequences of the test set samples, and the sequence similarity is calculated. The types and sources of mismatched bases are analyzed, and the accuracy of the reference database is evaluated. According to the evaluation results, the reference database is optimized and improved; for example, erroneous data are eliminated and new samples are added.

[0083] The detailed method for establishing the reference level database includes the following steps:

[0084] Step 1, select a DNA sample with a known sequence for nanopore sequencing to obtain raw current data;

[0085] Step 2: Process and analyze the current data to extract reliable level information;

[0086] Step 3, repeat 1-2 multiple times, accumulate level data, and establish a reference level database.

[0087] Step 4: Associate the levels in the reference level database with their corresponding known DNA sequences to form a "level-sequence" mapping relationship.

[0088] Step seven, based on the Needleman-Wunsch and Smith-Waterman sequence alignment algorithms, compare the level data obtained from the DNA sequencing to be tested in step five with the level data in the reference level database, find the most similar sequence, infer the DNA sequence of the DNA to be tested, add the level of the newly measured DNA sequence to the reference level database, and recalculate the median of different sequence combinations to update and improve the reference level database and improve its accuracy.

[0089] The alignment algorithm uses the Needleman-Wunsch and Smith-Waterman sequence alignment algorithms as the basis. The sequence measured in the database is A, and the real sequence measured is B. The Needleman-Wunsch or Smith-Waterman alignment of two base sequences A and B allows gaps in the two sequences. Due to possible gaps, the optimal alignment between the first nA base of A and the first nb base of B is one of the following:

[0090] The comparison scoring formula is as follows:

[0091]

[0092] in:

[0093]

[0094] Path optimization:

[0095]

[0096] 1.A's first Aseq i-1 All Bseq of bases and B j Optimal alignment between bases plus gaps in A;

[0097] 2. All Aseq of A i The first Bseq of base and B j-1 Optimal alignment between bases plus gaps in B;

[0098] 3. The first Aseq of A in two sequences i-1 The first Bseq of base and B j-1 The optimal alignment between the base plus the final matched base (within a window of 20) (such as formula 4) is implemented as follows: Fig. 9 The longer optimal alignment is recursively calculated based on the shorter optimal sub-alignment, and the entries in the alignment table are filled from the upper left corner to the lower right corner as shown. Fig.10 shown.

[0099] The sequence of levels sometimes goes backwards, and if the Needleman-Wunsch scheme is used directly, each entry in the alignment table will depend on the entry to the left due to the steps forward and on the entry to the right due to the steps backwards, which in turn will depend on the entry in question.

[0100] To address this, each step in the alignment is required to advance sequence A by one step. The alignment is therefore a best-case mapping of each level in sequence A to its corresponding level in sequence B, or to an empty level if no good match exists. The alignment is thus tracked past the skipped levels in sequence B, rather than skipping them.

[0101] A) The Needleman-Wunsch alignment considers both horizontal and vertical steps in the alignment with a penalty w corresponding to bases that do not match in one or the other sequence.

[0102] Winder steps, which indicate a matching base, usually have no penalty unless there is a mismatch.

[0103] B) Our nanopore arrangement forces each step forward along sequence A, but allows backtracking in sequence B. Affine penalties are assigned to jumps and backtracking: for example, backtracking by 3 levels will incur a backtracking penalty wback plus a backtracking penalty twice as large as wback again (as in Eq. 5).

[0104] C) Difference between the Needleman-Wunsch alignment and the Nanopore alignment when there is a jump in sequence B.

[0105] D) In ​​the following process, (tracing back through the global optimal mechanism, such as formula 6), we perform biological nanopore detection on the DNA sequence of the unknown species, obtain this level for database comparison, and after obtaining the final sequence, we can pass it into the shared database to annotate this species.

[0106] It should be understood that the above specific embodiments of the present invention are only used to illustrate or explain the principles of the present invention, and do not constitute a limitation of the present invention. Therefore, any modifications, equivalent substitutions, improvements, etc. made without departing from the spirit and scope of the present invention should be included in the protection scope of the present invention. In addition, the appended claims of the present invention are intended to cover all changes and modifications that fall within the scope and boundaries of the appended claims, or the equivalent forms of such scope and boundaries.

Claims

1. A data processing method for determining DNA nucleic acid organisms using a biological nanopore, further comprising the following steps: Step 1: According to the characteristics of the biological nanopore, the DNA sequence to be detected is divided into a combination of several bases, and then the DNA sequence to be detected is subjected to biological nanopore detection; and the original current change data of the DNA sequence to be detected is obtained; Step 2, correcting the invalid values ​​in the original current change data by using a correction algorithm; Step 3, performing high-frequency noise filtering on the corrected current change data; Step 4: divide the processed current change data into opening current and blocking current, and filter the events that meet the conditions; Step 5: perform level search and fitting for the events that meet the conditions. First, identify the level, determine the boundary of the current resistance level generated by different base combinations, and then remove the false level. Step 6: Establish a reference level database, select DNA samples with known sequences as test sets, perform nanopore sequencing and level extraction on them, compare the level data of the test set samples with the reference level database, obtain the annotated DNA sequence, compare the annotated DNA sequence with the known actual sequence of the test set samples, calculate the sequence similarity, analyze the type and source of the mismatched bases, evaluate the accuracy of the reference database, and optimize and improve the reference database according to the evaluation results, such as removing erroneous data and adding new samples; Step seven, based on the Needleman-Wunsch and Smith-Waterman sequence alignment algorithms, compare the level data obtained from the DNA sequencing to be tested in step five with the level data in the reference level database, find the most similar sequence, infer the DNA sequence of the DNA to be tested, add the level of the newly measured DNA sequence to the reference level database, and recalculate the median of different sequence combinations to update and improve the reference level database and improve its accuracy.

2. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 1, characterized in that: In the step 1, the DNA sequence to be detected is detected by biological nanopore using a sampling frequency of 500k.

3. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 1, characterized in that: In the step 2, the invalid value is a value that exceeds a preset percentage of the normal value, and the correction algorithm includes one of a median filter and a weighted average algorithm.

4. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 1, characterized in that: In step 3, a Bayesian low-pass filter is used to filter the high-frequency noise of the corrected original current change data. The sampling frequency of this project is 500KHz, the filter adopts Bessel 4th order, and the cut-off frequency is 50KHz.

5. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 1, characterized in that: In step 4, event refers to a blocking current event generated when DNA passes through the pore, and the screening of event includes the following steps: Step 1: First, based on a broad threshold detection, select the extracted event with a duration Time>1S and a current range between (0.1, 0.75) of the detected opening current; Step 2, initializing the local opening current value according to the number of detected events and the detected opening current; Step 3: Then analyze each local open state and each event one by one and further filter suitable events.

6. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 5, characterized in that: The event further screening includes the following steps: Step 1: Check the states before and after each numbered event in event: event_before and event_after, which represent the open states before and after the event respectively, and the initial value is NaN. Step 2: If event is not the first state, the median of the open state before event is calculated as event_before; if event is not the last state, the median of the state after event is calculated as event_after. Step 3: If the current event number is greater than 2, and there is data in the first two open states, try to use the exponential decay function to fit the data of the previous state and update it to event_before. Step 4: If the current event number is less than the open state number and there is data in the current state, try to use the exponential decay function to fit the data of the next state and update it to event_after. Step 5, exception handling, if event_before and event_after abnormally deviate from the 3sigma value of the open-state current detected by event, they are also set to NaN.

7. The data processing method for determining DNA nucleic acid organisms using a biological nanopore according to claim 1, characterized in that: The establishment of the reference level database in step 6 includes the following steps: Step 1, select a DNA sample with a known sequence for nanopore sequencing to obtain raw current data; Step 2: Process and analyze the current data to extract reliable level information; Step 3, repeat 1-2 multiple times, accumulate level data, and establish a reference level database. Step 4: Associate the levels in the reference level database with their corresponding known DNA sequences to form a "level-sequence" mapping relationship.

Images