DNA barcodes using a double nanopore system

The potential of a double nanopore system to determine DNA barcodes has been demonstrated experimentally. By carrying out Brownian dynamics simulation on a coarse-grained model DNA with protein tag (barcodes) at known locations along the chain backbone, we demonstrate that due to large variation of velocities of the chain segments between the tags, it is inevitable to under/overestimate the genetic lengths from the experimental current blockade and time of flight data. We demonstrate that it is the tension propagation along the chain’s backbone that governs the motion of the entire chain and is the key element to explain the non uniformity and disparate velocities of the tags and DNA monomers under translocation that introduce errors in measurement of the length segments between protein tags. Using simulation data we further demonstrate that it is important to consider the dynamics of the entire chain and suggest methods to accurately decipher barcodes. We introduce and validate an interpolation scheme using simulation data for a broad distribution of tag separations and suggest how to implement the scheme experimentally.

The use of digitized DNA-barcodes 1,2 in species identification [3][4][5] has been a standard technique for the preservation of Earth's biological diversities 6 . The extinction of species is not homogeneous across the globe, rather a strong function of location. Many species in tropical countries are declining rapidly being on the verge of extinction. The use of portable desktop equipment will be beneficial to carry out the tests locally in different countries bypassing the restrictions of bringing samples from one country to a laboratory located in another country. Nanopore based sequencing methods, such as, MinIon produced by Oxford nanopore 7 is an important step towards that goal which will eventually replace traditional Sanger's type of sequencing. Thus there is a genuine need to develop real-time on site desktop methods for in-situ but fast and accurate determination of genetic information contained in barcodes. Besides, there is a broad interest in applications of nanopore sensing beyond DNA-sequencing and barcode determinations, such as, discrimination, post-translational modifications, and unfolding of various proteins [8][9][10][11][12] .
A double nanopore platform Fig. 1a has demonstrated that it has the ability to outperform single nanopore based devices 13,14 in controlling the velocity and chain conformations during the translocation process. The captured DNA by both the pores is not only straightened by the tug-of-war forces present at each pore, but adjustable differential biases and feedback mechanism at each pore offers overall a better control on the translocation speed 15 . The most noteworthy aspect is the ability of multiple scans 13,16 of the translocating DNA through the pores by flipping the net differential bias that not only increases the statistical accuracy of measurements, but in principle capable of providing additional information about the physical processes, which to date are largely unknown, and requires a thorough theoretical investigation.
While a double nanopore system offers immense promises, preliminary experiments reveal that extracting genomic lengths can be complicated due to lack of experimental information about the entire chain. The experiments can extract information about the dwell time of the protein tags in each pore and the time of flight (TOF) from one pore to the other only (Fig. 1). These information would be sufficient to extract spacing between the barcode accurately if the entire chain was moving with the same velocity. As the protein tags experience larger frictional drags, it is expected that different parts of the chain will not move through the double nanopore system with uniform velocity which will then introduce error in genomic length determination. The dispersion in velocity will also depend on the magnitude of the differential bias forces ( f LR ), the pore width ( h pore ) and the distance between the two nanopores ( d LR ). Experimentally nanoscale technology poses challenges to vary these parameters and can often be expensive. However, these dependencies can be efficiently studied in a model system using computer simulations strategies [17][18][19] combined with theoretical methods developed in the last decade for studying translocating chain through nanopores [20][21][22][23] and validated by simulation studies 24,25 .
In this article, we present Brownian dynamics simulation results for a coarse-grained model of dsDNA with protein tags attached to it mimicking the essentials of the experimental setup. The simulation data shows that the velocities of the chain segments are indeed nonmonotonic as a function of the monomer index (Fig. 2). We further demonstrate that using the information about only the protein tags to extract barcode distances, as measured in an experiment, lead to an over/underestimation of the barcode distances. We then use the nonequilibrium tension propagation (TP) theory of Sakaue 23 to explain the non uniformity of the velocity profile. The underlying physical picture that emerges also provides clues for the under/over estimation of the barcode distances and direct us to an interpolation scheme to determine the barcode distances accurately.

Coarse-grained model and Brownian dynamics
Our coarse grained model consists of a polymer chain of 1024 beads with 8 protein tags translocating through a double nanopore system (Fig. 1a), inspired by the 48 kbp long double-stranded -DNA used in the experiment by Liu et al. 16 , where sequence-specific protein tags are introduced chemically so that the distance between any two tags (Fig. 1b) is known.
In the simulation, each bead (monomer) represents approximately 46 bp unit long dsDNA and one tag is roughly equivalent to 75 bp in the experiment which translates to approximately 2-3 beads. The mapping procedure is discussed in the Supplementary Information and the general scheme of the BD simulation strategy for a homopolymer translocating through a double nanopore has been discussed in our recent publication 17,18 . The protein tags are introduced by choosing the mass and friction coefficient at tag locations to be significantly different from the rest of the monomers along the chain (which requires modification of the BD algorithm as discussed in the Supplementary Information I and Supplementary Fig. S1). Instead of explicitly putting sidechains at the tag locations, we made the mass and the friction coefficient of the tags 3 times larger. Typically heavier biomolecules in a solvent are correlated with larger volumes which would produce larger drags. The volume effect is implicit in our scheme as this to a first approximation will renormalize the friction. This we find enough to resolve the distances among the tags. The location of the tags are also chosen in such a way so that www.nature.com/scientificreports/ genetic distances are disparate (some tags being close by and some are far apart) as shown in Fig. 1b and Table 1.
It is also noteworthy that there is no left to right symmetry so that the center of mass of the chain is not located at the center of the chain. The key question to answer is mimicking the double nanopore experimental protocol how accurately can one extract these genomic distances so that the method then can be applied to determine genetic lengths in unknown specimens.

Repeated scans and measurements
The measurement protocols are described in Fig. 1g-j. The differential bias (Fig. 1g). Once the last tag passes through the right pore, the bias is switched to 1h) so that the direction of translocation is reversed. Here, we report the case for | � Later we mention what happens when an asymmetry is present and | f LR | � = | f RL | . This process (called flossing-one flossing consists of one R → L and one L → R scan) is repeated 300 times and we record the experimentally measurable quantities-the dwell time and the time of flight (tof) as described next. We reserve the subscripts 1 and 2 for the left pore (pore 1) and right pore (pore 2) respectively. For sake of brevity we define quantities for the L → R only, implicating that DNA is translocation from left pore to right pore and the same definitions hold for the R → L translocation.
Dwell time. The co-captured dsDNA in a dual nanopore system provides two different ways of time measurements for the protein tags during translocation which can be translated to genomic lengths. Similar to a single nanopore setup, one can measure the dwell time in each pore for a double nanopore setup. For example, for the L → R translocation (Fig. 1i) The time of flight (tof). In addition to the dwell time measurements, in a double nanopore setup, one can also measure the time taken by a monomer as it leaves one pore and touches the upper boundary of the other pore during its voyage across the pore separation d LR (Fig. 1j). This is called the time of flight tof and defined as follows: Dwell velocity v D and time of flight velocity v tof . Accordingly, one can calculate both v D as well as v tof using Eqs. (1) and (2) as follows, www.nature.com/scientificreports/ Here �· · · � implies average over multiple scans through the left and the right pore which reduces the statistical error in the measurements.

Results
Non uniform velocity profile. Due to lack of L → R symmetry both v D and v tof have a positive slope along the direction of translocation due to propagation of the tension front 24 (Fig. 2a,b). Figure 2c shows that average of both the directions and the slope is close to zero as we are considering the symmetric differential bias | � f LR | = | � f RL | . As expected, the velocity distribution of the chain segments become non-monotonic with respect to the average velocity of the entire chain. The protein tags with heavier mass ( m tag > m bulk ) and larger solvent friction ( γ tag > γ bulk ) reside at the lower envelope of the graphs (Fig. 2c) while the dsDNA monomers reach their maximum velocity somewhere in between the tags. In addition to the average dwell and time of flight velocities v D and v tof , we have plotted the average velocity of the entire chain v chain , each represented as a solid line. For the choice of the parameters the dwell time velocities are larger and v D is about 20% larger than the v tof . This is a coincidence and not a generic feature. We have checked that keeping all the parameters the same, an increase in pore thickness h pore will result in an overall decrease of the dwell time velocities due to increased friction for a thicker pore. It is also worthwhile to note that although in simulation we can calculate the dwell time and time of flight for all the monomers and tags, experimentally these data (Eqs. 1 and 2) are measured for the protein tags only as the tags produce significant current blockades to be measured. However, the entire chain contributes to the dynamics of the tags, and it is this lack of information for the entire chain inevitably leads to inaccurate measurements of the genomic length as demonstrated below.
Barcodes from the segmental velocities. Let us first evaluate the consequence if we calculate the average velocity of the segment v L→R seg (m, n) connecting the tags T m and T n by approximating it to be the average velocity of the tags T m and T n only so that Here, (�τ ) L→R mn and (�τ ) R→L mn are the time difference of arrival (we call it tag-time-delay) of the m-th and the n-th tags at L/R pore during R/L → L/R translocation. Fig. 1k illustrates a specific case (�τ ) L→R 87 . Equation 6c provides the final distance d tof (m, n) averaged over multiple scans in each direction. Likewise, using dwell time velocity Eqs. (3a-3c), and replacing the subscript tof by D in Eqs. (5a,5b), one can derive equations analogous to Eqs. (6a-6c). A similar equation for the barcode distance using dwell time data from both the pores: The distribution of barcode distances with respect to T 5 using Eqs. (6c) and (7) are shown in Fig. 3a-d and summarized in the 4th and 3rd columns of Table 2 respectively. A closer look reveals the over/under estimation of the barcodes (columns 3 and 4) w.r.t the theoretical value (column 1) occurs when v D and v tof (columns 5 and 6) are greater or less than the average velocity of the entire chain v chain (column 7), and can be immediately discerned from Fig. 2. Furthermore this is an uncontrolled velocity approximation introduced in Eqs. (5a,5b) and depends on the contour length separation between the tags which is the unknown and yet to be determined. We further observe that since v tof (m) < v chain for all the tags m = 1, 8 , the barcode distances are underscored. On the contrary, v D (m) for the tags are more dispersed above and below v chain , and whenever v D (m) ≃ v chain , Eq. (7) gives a better agreement (for T 3 and T 4 ). If we replace the approximate velocities in Eq. (6c) by the constant velocity v chain of the entire chain this improves the estimates significantly. This is shown in column 8 (Barcode Method-I) and discussed later. We now explain the source of discrepancy using the non-equilibrium TP theory.
TP Theory explains the source of discrepancy and provides a solution. Unlike a rigid rod, tension propagation has a considerable effect in governing the motion of a semi-flexible chain in the presence of an external bias force. In TP theory 23 and its implementation in Brownian dynamics 24,25 the motion of the subchain in the cis side decouples into two domains. In the vicinity of the pore, the tension front affects the motion directly while the second domain remains unperturbed, beyond the reach of the TP front. www.nature.com/scientificreports/ In our case, after the tag T m translocates through the pore, following monomers are dragged into the pore quickly by the tension front, analogous to the uncoiling effect of a rope pulled from one end (Please refer to Supplementary Fig. S2 and the movie in the Supplementary Information). The onset of this sudden faster motion continues to grow (Supplementary Fig. S2a) and reaches its maximum until the tension front hits the subsequent tag T m±1 (for R/L → L/R translocation direction), having larger inertia and viscous drag. At this time (called the tension propagation time 25 ) the faster motion of the monomers begins to taper down to the velocity of the tag T m±1 . An example of this process is shown in Supplementary Fig. S2b. This process continues from one segment to the other and explains oscillatory characteristics in Fig. 2c. These contour lengths of faster moving segments in between two tags are accounted for neither in Eq. (6c) nor in Eq. (7). The experimental protocols are limited in extracting barcode information through Eqs. (6c) and (7) (measuring current blockade time) and could be the possible source of error.

How to estimate the barcodes accurately?
We now propose two methods that take into account the dynamics of the entire chain and correctly determine the barcodes and can be implemented in a dual nanopore setup experimentally. is the time delay of arrival at the pore between T 1 and T 8 for L → R translocation. We then estimate the barcode distance d L→R mn between tags T m and T n as The barcodes for the L → R and R → L translocation are shown in Fig. 3e. The average shown in Fig. 3f corresponds to column 8 of Table 2 as mentioned earlier. This is a significant improvement compared to the usage of the average tag velocity of the chain segments. This method will work if one can put two additional tags at known distances, at or close to the two ends of the DNA being scanned. Alternately, scan time information can be used to have a better estimate of the average velocity of the chain. In our simulation, we kept the scanning length L scan constant with starting and ending values ranging from 0.0976L to 0.902L. By using the constant scanning length L scan , the average scan velocity v scan can be used to determine barcode distances by replacing v chain ≃ v 18 in Eq. (9) with where τ scan (i) is the scan time for the ith event, and N scan = 300.
Method 2: Barcode using two-step method. Having gained a better understanding of the velocities of the monomers of the dsDNA segments in between the tags we now rectify Eqs. (6c) and (7) by taking a weighted average of the velocities of tags and DNA segment in between as follows. First, we estimate the approximate number of monomers N mn ≃ d L→R tof /�b l � ( b l is the bond length) by considering the tag velocities only using Eq. (6c). We further re-calculate the segment velocity accurately by incorporating weighted velocity contributions from both the tags and the monomers between the tags as follows: www.nature.com/scientificreports/ Here, n next are the number of neighboring monomers adjacent to the tags those share the same tag velocity. We checked that n next ≈ 1 − −3 does not make a noticeable difference in the final result. The barcode distances are finally calculated as for L → R translocation and repeating the procedure for R → L translocation, shown in Fig. 3g for both L → R and R → L translocation. The average shown in Fig. 3h corresponds to column 10 of Table 2. It is worth noting that (i) in Eq. (11) the tag velocities are more weighted and makes a difference when N mn is small, i.e., the contour length between the tags is small, in which cases the monomers in between the tags move with almost same velocity as that of the tags. In the other limit when N mn >> 1 , it is the chain velocity that dominates and one can safely ignore the velocity of the two tags (the 1st two terms in Eq. 11). Since the number of tags are only a few (8 in 1024 in our case), Eq. (9) works well excepting when the tags are close by. (ii) The "two step" weighting procedure in Eq. (11) is only approximate and has room for further improvement as one can interpolate from v tag to v chain with a suitable interpolation scheme.
Dwell time versus TOF. We have repeated the same protocol to correct the data from the dwell time measurement replacing v tof in Eq. (11) by v D listed in column 10 of Table 2. They are practically indistinguishable excepting for short tag distances.
Tag-time-delay matrix and the sum-rule. If we use v chain = constant to determine the barcodes as in Eq. (9), then the average tag-time-delay �(�τ ) mn � = 1 2 �(�τ ) L→R mn + (�τ ) R→L mn � will be proportional to the barcode distances. One can then form a heat map of the normalized tag time delay (�τ ) mn = (�τ ) mn /(�τ ) 18 in the form of a 8 × 8 matrix as shown in Fig. 4. The values in each square when multiplied by the appropriate scale factor ( (�τ ) 18 × v chain ) will reproduce the barcode distances and can serve as a nice visual about the relative distances between the tags. We find that indeed it reproduces barcodes of Fig. 1b as listed in Table 1.
Realistic pores and biases. Two pores in an experimental dual nanopores system are not exactly identical, about 5-10% differences in pore diameters are reported 16 . Likewise, when the differential bias is reversed for translocation directions L → R and R → L respectively. We ran additional simulations by offsetting the ideal conditions and checked that the same methods (Method I and II) work.

Summary and concluding remarks
Dual nanopore platform has immense promise and advantages compared to its single nanopore counterpart. In this Brownian dynamics simulation study, we mimicked an experimental platform and explained why extracting information from the tags only need to be corrected by taking into account the motion of the entire chain. We invoked TP theory to explain the nonuniform velocity distribution of the entire chain as a function of the monomer index. The protein tags introduce oscillation on the uniform velocity of the chain that depends on the tension propagation time from one tag to the other. We have checked that the information obtained from the time of flight data is more accurate compared to the dwell time data from the individual pores. We further discovered that the most reliable quantity is the tag-time-delay of the successive barcodes to arrive at the L/R pore. When the distance between the tags is large the tag-time-delay will straight translate to genomic length excepting for those cases when the tags are close by. Our two-step interpolation scheme will overcome this issue. This is due to the fact that roughly it is the average velocity of the entire chain and not the average velocity of the two tags that needs to be used to calculate the barcode distances. The heat-map of the normalized tag-time-delay provides and the corresponding sum-rule are the direct proof of the efficacy of this method. This study also indicates how to improve the measurement protocol. With some prior information about the tags if one can selectively attach heavier molecules at the tag positions so that they produce sharp dips on the velocity profile of the entire chain, then the procedure will be more accurate.
We would like to make some comments on relating these results with the time scale of the actual experiment and the effects of the pores. Our coarse-grained chain length will inevitably lead to a faster translocation process in simulation compared to the actual speed of the dsDNA. It is clear from watching the simulation movies that the tension propagation occurs in the model system well within the scan time. Thus, it is expected that in the actual experiment this will be more prominent due to the slower movement of the chain. Secondly, in experiments it is likely that there are cross-talks between the forces at each pore (the authors are indebted to Walter Reisner and William Dunbar for bringing this to their attention), while in our simulation the forces are strictly local. Moreover, experimental pores are much wider, and thus folded configurations can translocate which are not included in our model with narrow pores. The tags can introduce additional pore-poymer interaction. These items will not change the main conclusion of the paper, however, will be included in our subsequent studies of polymer translocatoin through a double nanopore system in a broader context. We believe that this study will be immensely useful for designing future double nanopore platforms so that the data in the time domain can be translated to unravel fine structures of genomic lengths.