- Email: [email protected]
A mathematical model of recombinase polymerase amplification under continuously stirred conditions
Biochemical Engineering Journal 112 (2016) 193–201
Contents lists available at ScienceDirect
Biochemical Engineering Journal journal homepage: www.elsevier.com/locate/bej
FA mathematical model of recombinase polymerase amplification under continuously stirred conditions Clint Moody a , Heather Newell a , Hendrik Viljoen a,b,∗ a b
Department of Chemical and Biomolecular Engineering, University of Nebraska-Lincoln, Lincoln, NE 68588, USA Preclinical Drug Development Platform, North-west University, Potchefstroom 2520, South Africa
a r t i c l e
i n f o
Article history: Received 18 January 2016 Received in revised form 18 March 2016 Accepted 19 April 2016 Available online 23 April 2016 Keywords: Recombinase Polymerase Amplification (RPA) Isothermal PCR Couette flow Mathematical model Disease diagnostics
a b s t r a c t Growing interest surrounds isothermal PCR techniques which have great potential for miniaturization for mobile diagnostics. Particularly promising, Recombinase Polymerase Amplification (RPA), combines this advantage of isothermal PCR with simplicity and rapid amplification. A mathematical model is presented of Recombinase Polymerase Amplification (RPA) and compared to experimental data. This model identifies the rate limiting steps in the chemical process, the effects of stirring, and insights in to using RPA for quantitative measurement of initial DNA concentration. Experiments are shown in which DNA amplification occurs under conditions of Couette flow and conditions of rotational turbulent flow. Hand mixing has been shown to dramatically shorten amplification times but introduces unpredictable variability. In some cases, this variability manifests itself as human error induced false negatives, a serious problem for all potential applications. Mechanical stirring demonstrates similarly short delay times while retaining high repeatability and reduces the potential for human error. Published by Elsevier B.V.
1. Introduction Polymerase Chain Reaction (PCR) has revolutionized molecular biology and has proven itself as one of the fastest and most specific methods of disease detection. In many places around the world, the sparsity of hospitals and clinics make developing point of care procedures a necessity. For PCR, the requirement for bulky and expensive thermocyclers remains one of the main obstacles for point of care use. To address this issue, many are turning towards isothermal PCR for the answer. A technique that does not require thermocycling to denature double stranded DNA but instead only requires a thermal bath lends itself much more readily to point of care diagnostics. Several isothermal techniques have emerged since the development of PCR. The most widely used are rolling-circle amplification, helicase-dependent amplification (HDA), Loop Mediated Isothermal Amplification (LAMP), and Recombinase Polymerase Amplification (RPA). Amplification times for Rolling circle generally take 65 min [1]. In addition to this 65 min process, it also requires significant preparation steps beforehand that can significantly drive up the total time
∗ Corresponding author. E-mail addresses: [email protected] (C. Moody), [email protected] (H. Viljoen). http://dx.doi.org/10.1016/j.bej.2016.04.017 1369-703X/Published by Elsevier B.V.
and overall complexity. The simplest and most efficient method for generating a suitable circular template, using a padlock probe, requires carefully designed oligomers and limits the amplification to small target regions [2]. Helicase-dependent amplification employs a helicase to unzip the double stranded DNA template to accommodate primer annealing. Helicase activity, however, requires bubbles to form naturally at AT rich regions [3]. The high GC content in Mycobacterium and other bacteria, in many cases, would inhibit this technique. Another drawback of HDA is that optimal amplification only occurs for short sequences of 80–120 bp [4]. Loop Mediated Isothermal Amplification (LAMP) has the benefit of being a simple method to perform from a kit but is deceptively complicated. It requires 4 primers that need to bind at 6 locations. The DNA regions around the primers must form loops without hairpins or other secondary structure so primer design is no trivial task [5]. Even though it does not require a high rate of thermocycling, it does require an initial denaturation step at 95◦ Celsius before amplifying at 65◦ C [6]. Amplifications with this method generally require 30–60 min [6]. A newly developed method called Recombinase Polymerase Amplification (RPA) shows great promise especially for pointof-care diagnostics. RPA uses recombinase to insert primers in to double stranded DNA rather than denaturing DNA using high temperature cycling [7]. This allows RPA to operate at a con-
194
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Nomenclature G m R FG n H3 PO4 FR FR FnR FnR D FD P PFD B np Sf t1 2
k1 k1f
Single stranded DNA binding protein (Gp32) Number of Gp32 binding sites on a primer Recombinase complex (UvsX.6*UvsY) Forward primer/Gp32 complex Number of recombinase binding sites on a primer Inorganic phosphate Unstable forward primer/recombinase complex Stable forward primer/recombinase complex Forward primer complexed with unstable filament of n recombinase molecules Forward primer complexed with unstable filament of n recombinase molecules DNA template Forward primer/DNA complex Polymerase Polymerase/forward primer/DNA complex Number of base pairs in template Number of base pairs in primer Empirical scaling factor to modify diffusion limited rate constants Time required for DNA concentration to reach 50% of saturation Equilibrium constant for Eq. (1) 105 M −1 Forward rate constant for Eq.
108 1/ (M × s)
5×
(1)
k3f
Reverse rate constant for Eq. (1) 5 × 103 1/s Equilibrium constant for Eq. (2a) 68 × M −1 Forward rate constant for Eq. (2a) 108 1/ (M × s) Reverse rate constant for Eq. (2a) 1.471/s Rate constant for Eq. (2b) 47 × 10−3 1/s Equilibrium constant for Eq. (2c) 3 × 106 M −1 Forward rate constant for Eq. (2c) 108 1/ (M × s) Reverse rate constant for Eq. (2c) 33 1/s Rate constant for Eq. (2d) 4.6 × 10−3 1/s The Michaelis constant for Reactions (3) and (4a) 20.35 × 10−6 M Forward rate constant for Eq. (3) 108 1/ (M × s)
k3r
Reverse rate constant for Eq. (3) 59.37 1/s
k1r keq2a keq2af keq2ar k2b keq2c keq2cf keq2cr k2d KM3a
k4acat k4bcat k5f
Forward rate constant for Eq. (4a) 4.22 1/s Forward rate constant for Eq. (4b) 8.32 1/s Forward rate constant for Eq. (5)
107 1/ (M × s)
2. Materials and methods 2.1. Amplicon design for real-time results 1.2 ×
k5r
Reverse rate constant for Eq. (5) 0.06 1/s
k6f
Forward rate constant for Eq. (6) 87 1/s
k7
Forward rate constant for Eq. (7) 4.1 1/ (M × s)
k8
Forward rate constant for Eq. (8) 1.13 1/ (M × s)
time PCR [7]. While many detailed models have been presented for the other isothermal PCR methods mentioned, RPA differs from them in some key ways that warrant special attention. To compare the model with data, the experimental methodology must minimize the effect of poorly repeatable inputs. One example of this would be hand mixing. The current methodology recommended by the supplier of RPA technology, Twistamp, requires two hand mixing steps. Our experience has shown that these mixing steps introduce a significant amount of variability in results. The tubes must be mixed both at the start of the process and after 4 min. The removal of the tube for mixing at 4 min creates a gap in data collection at a crucial point in the amplification process partly defeating the purpose of real time detection. To make matters worse, if the tubes are insufficiently mixed at 4 min, the samples will not amplify or amplification will be stunted. Another issue arises when, during the mixing step, droplets stick to the lids of the tubes artificially reducing the signal levels. These sources of human error make having a skilled operator imperative and confounds efforts to compare results obtained by different investigators. While these mixing steps are error prone and introduce significant randomness, without them the process takes several times as long to give a result. In this paper, we will demonstrate an automated stirring process that eliminates these sources of human error and simplifies the process to make it more easily employed in field applications. Automated stirring has been demonstrated before using lab-on-a-chip technology [10], but our method is more accessible to a wide range of users without buying specialty disposable instrumentation. Initially, our procedure used an ordinary benchtop drill press to stir the reaction mixture. The equipment was then scaled down to use instead a small RC motor. The only part of the set-up that needed replacement after each test was the drill bit/stir rod. Not only does this motor take up negligible space (about as much as a 1.5 ml tube) but its low cost also makes an array of 8, one for each tube, quite feasible even in low budget applications. This stirring method, by removing the arbitrary element of hand shaking and moving the reaction in to a kinetically controlled regime, made the reaction more conducive to mathematical modeling. We herein present, to our knowledge, the first mathematical model specifically addressing the unique considerations of RPA.
stant temperature of around degrees celsius—much lower than other isothermal techniques like LAMP [8]. Using the tubescanner from TwistDx, the instrumentation costs for RPA are significantly cheaper than those for either PCR or LAMP and include fluorometers for real time detection [8]. It only requires two primers that, with a length of 30–35 base pairs, are only modestly larger than PCR primers [9]. RPA also boasts shorter amplification times than LAMP. In combination with real time fluorescent measurements, a result can be determined in approximately 10 min [7]. Remarkably, the rapid rate of amplification does not come at the expense of reduced sensitivity. RPA can amplify as few as 10 initial copies of template—a detection threshold similar to that of traditional real
DNA was extracted from Mycobacterium smegmatis and Escherichia coli (as a negative control) using the phenol-chloroform method. The amplicon for the detection of the Mycobacterium Smegmatis DNA was designed using the sequence IS1096. The RPA amplicon was 116 base pairs long. SYBR green was used instead of designing a specific probe. The forward and reverse primers were: IS1096F: 5 -CTCATCGAACATTCCCGCGAACACGTTCCGACCAG-3 IS1096R: 5 -CTTGACGGTGTAGAGACGATCAGCTGCTTTCGC-3 2.2. RPA conditions and detection RPA was performed using the TwistAmp Basic kit, having a 50 l volume. Reaction mixtures were formulated according to TwistDx Ltd. recommendations with one of notable change. Instead of using probes for detection, 5× concentrated SYBR Green was added to each reaction mixture. The DNA concentration was standardized at 0.1 ng/l for all tests. Once all components were added, the tubes were vortexed and centrifuged; a solution of 14 mM Mg acetate was then pipetted into each reaction tube which was then immediately put into the tubescanner device. The reaction components were
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
supplied by TwistDx Ltd., Cambridge, United Kingdom. The fluorescent signal is in the FAM channel, 470 nm excitation and 520 detection. Each experiment was run for 20 min at 39◦ C. The realtime readout was courtesy of Twista Studio software, version 2.6.0 (TwistDx, Ltd.) which interpreted fluorescent data into graphs of Intensity (mV) vs. time (s). 2.3. Experimental set-up For all reactions stirred at or below 3050RPM, RPA reaction mixtures were stirred with a PerformaX 10” benchtop drill press using 1/16th in. stainless steel drill bits. Following each experiment, the drill bit was replaced to avoid contamination. Experiments using 1/16th in. diameter polypropylene rods (Data not shown) gave similar results but required higher rotation speeds. In order to allow the drill access to the reaction mixture, the lid of the TwistA system had to be removed. Since this exposed the optical detection system to ambient lighting, all experiments were performed with the lights out. For our highest stirring speed test, at 7050RPM, we used an RC motor. The exact model of motor was the Ares AZSH1313 Ares Ethos QX130 Motor w/Pinion Gear and Wire Leads, Clockwise Rotation. A brass coupler and two set screws were used to mount the drill bits to the RC motor. This particular motor becomes unreliable at lower speeds and thus could not be used in the tests below 7050RPM. In many applications, high speed stirring provides the best results. Only when used with oligomer probes do we see improved performance with slower stirring speeds. For these applications, other similar motors with slower optimal speeds could be substituted. 3. Theory PCR employs thermal cycling to denature double stranded DNA and facilitate primer annealing. Recombinase polymerase amplification, instead, uses a recombinase complex to denature double stranded DNA so that primers can bind. This recombinase complex consists of an UvsY hexamer bound to UvsX. In our model, we assume that UvsY forms a stable hexamer and binds to UvsX prior to the start of the test. Under normal RPA conditions, this is a safe assumption because the protein pellet contains both proteins and the UvsY hexamer is stable with a small dissociation constant [11]. Another essential element of the RPA pellet, a single stranded DNA binding protein called Gp32, helps remove secondary structure from primers to facilitate the binding of recombinase [12]. We assume that Gp32 binding must proceed recombinase binding. Thus the first reaction treated by our model involves primers binding reversibly with Gp32, represented simply by G. m × G + F FG
(1)
For the sake of simplicity, we assume perfect symmetry between the reactions involving forward and reverse primers. This allows us to reduce the number of differential equations in the model by simply multiply the appropriate terms by a factor of two. The RPA pellet composition as a patented technology is publically available and, assuming a 50 ul reaction, we calculate a 9.4 × 10−7 Mconcentration of Gp32 [13]. The published value for the equilibrium constant of this reaction, k1 , is 105 M −1 [14]. The dissociation rate k1r is5 × 103 1/s [14]. The binding site for Gp32 is 7 base pairs long so for a primer 32 base pairs long m is roughly 4 [14]. Next, the recombinase complex, R, displaces Gp32 from the primer in a two-step process [15]. First, the process requires a nucleation step. R + FG FGR
FGR → FGR + G
(2a) (2b)
195
The complex binds to the primer to become the unstable complex FGR but only becomes stable when it displaces Gp32 to form FGR [15]. While UvsX uses ATP hydrolysis to drive this process, we assume that ATP levels are constant and not limiting. The presence of phosphocreatine kinase and phosphocreatine in the pellet acting to regenerate ATP makes this assumption justified. Once one UvsX.UvsY complex has bound and achieved stability, nucleation is complete and the process can proceed to the growth phase. FGR + (N − 1) × FGnR
(2c)
FGnR → FnR(m − 1) × G
(2d)
Each UvsX monomer binds to 4 base pairs of DNA so the number of monomers consumed, n, is approximately 8 for primers of length, np = 32 [16]. Liu et al. lists the kinetic parameters and equilibrium constants for these four reactions above [15]. This paper presents alternative values using either hyperbolic or global methods. The global parameters better fit the data so they were used in all four reactions. The value of the equilibrium constant of Reaction (2a), keq2a ,as determined by global methods is 68 × M −1 [15]. The value of the forward rate constant of reaction (2b), k2b , is 47 × 10−3 1/s [15]. The value of keq2c , another equilibrium constant this time for Reaction (2c), for kinetically controlled reactions is 3 × 106 M −1 [15]. The rate constant for the complete dissociation of Gp32, k2d , is 4.6 × 10−3 1/s [15]. Based on the pellet composition the concentration of UvsX is 5.9 × 10−6 M [13]. Once recombinase has bound a primer, it seeks out homology in the DNA template. It then inserts the primer in to the DNA and dissociates from both by hydrolyzing ATP. We originally modeled this reaction using Michaelis-Menten kinetics as in Farb et al. We quickly discovered, however, that the underlying assumptions of the Michaelis-Menten kinetics poorly described this system. So we instead used the more fundamental two reaction formulation. FnR + D + 2 × m × G FnRD
(3)
FnRD + n × ATP → FD + n × R + n × AMP + n × PPi
(4a)
In Reaction (3), the recombinase filament escorts the primer to the region of homology in the DNA template, here represented by D. Next, hydrolysis of ATP drives the dissociation of the recombinase from the primer/DNA complex, FD. Michaelis-Menten kinetics assumes that the equilibrium reaction has reached steady state. This assumption in particular does not apply during the transient when diffusion rates are slow and substrate concentrations are low. We used the constants k4acat = 4.22 1/s and, for the unstirred reaction, KM3a = 20.35 × 10−6 M [16]. Along with the assumption that the forward reaction rate k3f = 108 1/ (M × s), corresponding to the diffusion limited case, we were able to derive the rate constants from the Michaelis-Menten parameters. This gives us a reverse reaction rate, k3r, of 59.37 (1/s). The binding of GP32 to the opposing DNA strand helps stabilize the complex formed in Reaction (3). Unlike RecA and other members of the recA family, UvsX has two modes of ATP hydrolysis. It can hydrolyze ATP to either AMP or ADP [16]. As side products, Reaction (4a) produces pyrophosphate (PPi) and Reaction (4b) produces inorganic phosphate (H3 PO4 ).To account for this secondary mode of ATP hydrolysis we added the following reaction. FnRD + n × ATP → FD + n × ADP + H 3 PO4
(4b)
We used the forward rate constant k4bcat = 8.32 1/s for this reaction [16]. RPA employs strand displacing Sau Polymerase (a polymerase derived from Staphylococcus aureus denoted here as P) to extend primers [7]. Polymerase can bind reversibly to DNA primer complexes. We assume binding rates comparable to those
196
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 1. Overall reaction scheme diagram of RPA.
of DNA polymerase I which are well established as k5f = 1.2 ×
107 1/ (M × s) and k5r = 0.06 1/s [17]. FD + P PFD
(5) −6
The concentration of Polymerase is 1.34 × 10 M [13]. Once the tertiary complex has been formed polymerase can extend the primer to form a new strand of DNA using nucleotides (dNTP’s). The factors of two in this equation indicate that a simultaneous extension of the reverse primer generates a second strand of DNA. Again we must assume a rate constant of similar magnitude as E. coli Polymerase I which has been measured to be approximately k6f = 87 1/s [18]. PFD + 2 × B × dNTP → P + 2 × B × PPi + 2 × DNA
(6)
If pyrophosphate exists in sufficiently high quantities, it can drive the reverse reaction [19]. The reverse reaction rate was calculated from the forward reaction rate and the Gibbs free energy change of the extension process. Using the Gibbs free energy reported in Kuchta et al., the calculated the equilibrium constant is approximately 77 [17]. Thus k8 = k6f /77. Since this degradation is proportional to the length of DNA k7 = k8 × B/32. This process had little effect on the model outcome under the experimental conditions examined but may be significant for some reactions. FD + np × PPi → D + np × dNTP
(7)
D + B × PPi → B × dNTP + FD
(8)
This paper focuses mainly on the TwistDx Basic and FPG kits. These pellet formulations do not include exonuclease. The equations, however, can be easily modified to accommodate its presence. Reactions for exonuclease activity would take the same form as (7) and (8) (Fig. 1). Before presenting the rate equations for RPA, it is important to note that the pieces of the model came from multiple papers and we needed to make the different formulations compatible. We define DNA concentration on a per copy basis but our sources define it based on the number of binding sites or the number of base pairs. The rate constants for the Gp32 equilibrium used Gp32 binding sites as their reference. The rate and equilibrium constants for the recombinase reactions were calculated per recombinase binding site. Thus we need correction factors for our per primer calculations. The rate constants for primer insertion in to dsDNA were defined based on a per base pair concentration of DNA rather than the copy number so correction factors account for this difference in methodology (Fig. 2). The rate equations we use are: d [R] = −2 × keq2af × [R] × (n2 × [FG]) + 2 × keq2ar dt
× FGR − 2 × keq2cf × [R] × FGR + 2 × keq2cr × FGnR
+ 2 × k4acat × [FnRD] + 2 × k4bcat ∗ [FnRD]
Fig. 2. Model predictions for DNA concentration compared with three replicates of experimental data with a stirring speed of 7050RPM.
d ([FG]) = 1/m × k1f × [G] × [F] − 1/m × k1r × [FG] − keq2af dt
× [R] × (n × [FG]) + 1/n × keq2ar × FGR )
d [FGR ] = keq2af × [R] × (n × [FG]) − 1/n × keq2ar dt
× FGR − k2b × FGR
d [FGR] 1 = k2b × FGR − × keq2cf × [R] × [FGR] n−1 dt +
1 × keq2cr × FGnR n−1
1 d [FGnR ] × keq2cf × [R] × [FGR] = n−1 dt −
1 n−1
× keq2cr × FGnR − k2d × FGnR
d [FnR] = k2d × FGnR + dt ×k3f × [FnR] × [D]
1 n
× k3r × [FnRD] −
1 n
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 3. Model predictions for the concentration of primers forming an unstable complex with one molecule of recombinase. d[FnRD] dt
=
1 n
k4acat × [FnRD] − d[D] dt
=−
1 n
× k3f × [FnR] × [D] −
1 n
1 n
× k3r × [FnRD] −
n
Fig. 4. Model predictions for the concentration of primers forming a stable complex with one molecule of recombinase.
×
× k4bcat × [FnRD]
× k3f × [FnR] × [D] +
1 n
× k3r × [FnRD] + 2 ×
k6f1 × [PFD] + 2 × k7 × [PPi] 1 × [FD] − 2 × k8 × [PPi] ∗ [D] n
1
197
× k4acat × [FnRD] +
n
d[FD] dt
=
× k4bcat × [FnRD] − k5f × [P] ×
[FD] + k5r × [PFD] − k7 × [PPi] × [FD] + k8 × [PPi] × [D] −k5f × [P] × [FD] + k5r × [PFD] + k6f × [PFD]
d[P] dt
=
d [PFD] = k5f × [P] × [FD] − k5r × [PFD] − k6f × [PFD] dt d [dNTPs] = −2 × B × k6f × [PFD] + 2 × B × k8 × [PPi] dt × [D] + 2 × np × k7 × [PPi] × [FD]
d [PPi] = 2 × k4acat × [FnRD] + 2 × B × k6f × [PFD] − 2 × B dt ×k8 × [PPi] × [D] − 2 × np × k7 × [PPi] × [FD]
4. Results and discussion Using the model described above, we compared with experimental results and found reasonable agreement with data taken at 7050RPM. This best corresponds to one of the fundamental assumptions of this model: namely that the fluid is well mixed and all reactions are kinetically controlled. The Reynold’s number of this rotational flow is roughly 5000 so turbulence will dramatically enhance the mixing of the reactants. Without consistent mixing to minimize diffusion effects, there would be no foundation for this model. Later, we will relax the kinetically controlled assumption to describe the effect of slower stirring rates but first we must characterize the simplest form of the model. To understand the time evolution of the RPA reaction, we simulated the differential equations using a second order Euler method. Our simulations showed that the primers and Gp32 bind almost completely and instantaneously. So, for appropriate Gp32 concentrations, one could safely assume that all primers exist as this complex at t = 0s with negligible error. Likewise, FGR’ appears almost instantaneously (Fig. 3).
Fig. 5. Model predictions for the concentration of primers forming an unstable complex with a filament of n molecules of recombinase. In this case n = 8.
The first step to take any significant time (2b) produces peak levels of FGR at around 50s (Fig. 4). The production of FGnR’ (2c) generally tracks the concentration of FGR but introduces significant smoothing (Fig. 5). Fig. 6 shows that the rate limiting reaction for RPA is Reaction (2d) which suggests that the process could be improved by reducing the binding affinity of the single strand binding protein. The time evolution of the next three species produced (FnRD, FD, and PFD) appear qualitatively similar to the concentration of FnRD in Fig. 7. Under the conditions of our experiments, not only are they qualitatively similar, they roughly proportional. Estimating PFD as 0.475*FnRD, only results in a 4% error in FnRD across the peak in concentration (Fig. 8). Thus we could simplify the model by eliminating two differential equations while generating nearly identical results. This estimate is valid for a wide range of initial DNA concentrations so long as other variables such as polymerase concentration remain constant. For polymerase concentrations below optimized levels, FnRD and FD remain roughly proportional but limited polymerase concentrations control the production of PFD.
198
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 6. Model predictions for the concentration of primers forming a stable complex with a filament of n molecules of recombinase. In this case n = 8.
Fig. 8. Model predictions showing that the concentration of polymerase/primer/DNA complex is roughly proportional to the concentration of primers being delivered to the DNA template.
Fig. 7. Model predictions for the concentration of recombinase delivering primers to the DNA template.
The key reaction in the process is (3). The limiting reactant in this process controls the rate of DNA production. RPA, for the conditions investigated, occurs in three distinct stages. During the first stage, the assembly of the required FnR occurs. Once these primer/recombinase complexes form and far exceed the number of available DNA molecules, the second stage begins. The second stage shares much in common with the various methods of isothermal PCR. During this stage, the limited supply of DNA determines the rate and DNA concentration follows a roughly exponential curve. Once DNA reaches high enough levels, FnR once again becomes rate limiting. Without measurements of FnR, it is impossible to see the transition between the first and second stages based on a DNA plot alone (except in the unusual circumstance that DNA concentrations begin at detectable levels). The transition between the second and third stages, however, can be seen as a rounded corner in the graph which ends the exponential growth phase. These three stages result in behaviour that differs qualitatively from traditional PCR. Since FnR production occurs independently
Fig. 9. Model predictions of the effect of varying polymerase concentrations on the rate of DNA production.
from reactions after (3), any changes after that point will only affect phase 2. If some quantity such as polymerase concentration changes, the graphs continue to converge at the beginning of stage 3. Changes made that effect reactions before (3), however, affect all three stages. The major contribution to the delay time, under our experimental conditions, is the extension step of the recombinase filament assembly. This contrasts sharply with the time evolution in other forms of isothermal PCR. While the graphs may look qualitatively similar, the delay time in RPA is much less sensitive to the initial DNA concentration. As can be seen from Fig. 9, the polymerase concentration has been optimized for the commercial pellets and any further increases in concentration or improvements in the polymerase processivity would have almost no effect. Decreasing polymerase concentration, on the other hand, can significantly increase the
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 10. Model predictions of the effect of different initial DNA concentrations on DNA production.
199
Fig. 11. Model predictions showing increased sensitivity to DNA concentration when Polymerase concentration is reduced by a factor of 100.
delay time. This second effect may find unexpected usefulness in quantitative studies of initial DNA concentration. As has been shown previously, using quantitative PCR techniques one can arrive at a reasonable estimate of the initial concentration of DNA simply extrapolating backwards assuming exponential growth. This same method naively applied to RPA will result in nonsensical estimations orders of magnitude smaller that the true initial value. For the current optimized pellet we must instead use the shifted equation: t 1 = −14.5 × log10 (D) + 16.5 2
where t 1/2 is the time required for the signal to reach ½ its maximum value. The time shift accounts for the time required to assemble the recombinase filament on primers. Fig. 10 shows the effect of different initial DNA concentrations and one immediately notices that in its current form RPA is not terribly conducive to quantitative studies because such small differences in delay time can easily fall below the experimental variability. Reducing polymerase concentration or employing a less processive polymerase could potentially resolve this problem. Both methods would slow DNA extension while leaving the RPA specific reactions relatively unaffected. Since commercial pellets contain both the polymerase and the recombinase, it may be advisable to reformulate the pellet composition for quantitative RPA. Simply reducing the quantity of the pellet would reduce the recombinase levels below the concentration needed for the reaction to proceed. To demonstrate this spreading effect DNA dilutions were simulated at the current Polymerase concentration as well as at a hundredth of that concentration (Fig. 11). For the optimized polymerase concentration the delay time to half maximum is given by the formula above. For the reduced polymerase case, the formula is instead:
Fig. 12. Model predictions for the effect of varying recombinase concentrations on the rate of DNA amplification.
dramatic effect on the overall amplification. Amplification drops to nearly zero at 1/3 of the required recombinase. At these low concentrations, primers must compete for scarce recombinase complexes and complete filament assembly becomes exceedingly rare. Thus the entire process becomes stalled between the nucleation step (2b) and the growth step (2c). 4.1. Experimental stirring
t 1 = −33.2 × log10 (D) + 98.4 2
The spreading increased by a factor of 2.3 while the delay has only increased by a factor of 1.7. Using the right pellet formulation, this method could potentially generate quantitative results of any desired resolution. The effect of varying concentrations of recombinase was also simulated. From Fig. 12, we can see the insensitivity of DNA amplification to excess recombinase. Once recombinase levels become insufficient, however, small differences in concentration have a
Initially, a concern with mechanical stirring of the reaction mix was that the drill bit would interfere with the optical signal. Fortunately, experiments showed that a 1/16th in. drill bit submerged 3.5 mm in to the fluid had only a small effect on the signal produced by a 50 l reaction volume. For larger diameter drill bits, a constant offset was observed which could easily be subtracted from the signal to reproduce the expected results. While TwistDx recommends using specialty DNA oligomer probes, we did not observe any issues with using SYBR green either in terms of amplification
200
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
time for all stirring speeds investigated. The equation of the best fit line with an R2 value of 0.51 was: t 1 = −0.0082 × RPM + 273 2
Fig. 13. Experimental results for different mixing rates averaged over 3 replicates.
Tests done with stirring speeds below 570RPM amplified but their time evolution displayed poor repeatability (Data not shown). This poor repeatability was especially pronounced in experiments conducted with no mixing at all. While the rate of reaction can be dramatically improved by hand mixing, this introduces significant problems for repeatability as well. Aside from the inherently arbitrary nature of hand mixing, care must also be taken so that droplets do not stick to the lid of the tubes. These droplets contain DNA that will be hidden from the detector and thus can greatly diminish the observed fluorescence. Most significantly, we observed a higher rate of false negatives with hand mixing than with drill stirring. For any application, false negatives should be kept to a minimum and any attempt at quantifying the initial DNA concentration using RPA will require a mixing method with highly repeatable outcomes. Now we are ready to discuss modifications to the model that incorporate different stirring speeds. The forward reaction rates of (2a), (2c), and (3) are all diffusion limited. Thus we expect the corresponding reaction rate constants to decrease as stirring speed decreases. To model this we introduce a scaling factor, Sf, to reduce the forward rate constant for each of these reactions. Based on our experimental data, the appropriate Sf can be roughly calculated as:
Sf = 8.05 × 10−5 × RPM + 0.347 To use this scaling factor, simply replace keq2af, keq2cf, and k3f with Sf*keq2af, Sf*keqcf, and Sf*k3f respectively. We determined the scaling factor for each experimental condition by generating a simulation with the same delay time as our data collected at the different stirring speeds. This represents merely a simple empirical formulation which may prove useful to the experimentalist. The underlying chemical phenomenon cannot be determined from our limited data but clearly involves the enhancement of diffusion controlled Reactions—(2a), (2c), and (3) being the most significant. Since this data was taken at Reynolds numbers ranging from roughly 400–2200 (all within the regime of laminar flow), it should not be extrapolated above 3200RPM at which point the flow begins to transition to turbulent flow. Fig. 14. Experimental results showing the relationship between the time required for the reaction to reach half its maximum value and the stirring speed in RPM.
5. Conclusions inhibition or increased noise. Since probes must be tailored specifically for each application which represents a significant cost, and must be optimized empirically, we recommend using SYBR green in applications that do not require multiplexing and exhibit minimal nonspecific amplification. Using E. coli DNA, we showed that no detectable non-specific amplification occurred in one hour of incubation. So, unless otherwise noted, all experiments were conducted with SYBR Green however similar results were obtained with the recommended oligomer probes. The only significant difference we observed between the probes and SYBR Green was that when drill speeds exceeded 900RPM it compromised the probes’ ability to fluoresce. Probes appeared to work normally at drill speeds of 570RPM so the system will work with probes given the proper choice of drill speed. Increasing rotation speed increases the rate of amplification (Fig. 13). The relationship between stirring speed measured in RPM and the delay time of the amplification is roughly linear (Fig. 14). We continue to see a small but consistent improvement in delay
1. The equations listed represent to our knowledge the first mathematical model presented for Recombinase Polymerase Amplification 2. Consistent quantitative experimental results can be achieved using mechanical stirring eliminating human variability as well as making the process more reliable. 3. Stirring significantly increases the rate of amplification. 4. Insufficient recombinase dramatically hinders or even stalls the reaction. 5. Amplification can be understood as occurring in 3 stages: the first and third limited by FnR while only the second is limited by DNA concentration. The presence of the first stage makes estimating the initial DNA concentration different than the methods used in traditional quantitative PCR. Exponential growth models can underestimate initial DNA concentration by several orders of magnitude. 6. The model gives insights of value to anyone attempting to use RPA to quantify DNA concentration.
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Acknowledgements We thank D.K. Zinniel for providing DNA samples from Mycobacterium Smegmatis and E. coli. References [1] Chia-Chen Chang, et al., Diagnostic devices for isothermal nucleic acid amplification, Sensors 12 (6) (2012) 8319–8337. [2] Jia Li, Joanne Macdonald, Advances in isothermal amplification: novel strategies inspired by biological processes, Biosens. Bioelectron. 64 (2015) 196–211. [3] Schaerli Yolanda, et al., Isothermal DNA amplification using the T4 replisome: circular nicking endonuclease-dependent amplification and primase-based whole-genome amplification, Nucleic Acids Res. 38 (22) (2010) e201. [4] Kong Huimin, Myriam Vincent, Yan Xu, Helicase dependent amplification of nucleic acids. U.S. Patent No. 7, 282, 328. 16 Oct. (2007). [5] Notomi Tsugunori, et al., Loop-mediated isothermal amplification of DNA, Nucleic Acids Res. 28 (12) (2000) e63. [6] Tomita Norihiro, et al., Loop-mediated isothermal amplification (LAMP) of gene sequences and simple visual detection of products, Nat. Protoc. 3 (5) (2008) 877–882. [7] Euler Milena, et al., Recombinase polymerase amplification assay for rapid detection of Rift Valley fever virus, J. Clin. Virol. 54 (4) (2012) 308–312. [8] Teoh Boon-Teong, et al., Early detection of dengue virus by use of reverse transcription-recombinase polymerase amplification, J. Clin. Microbiol. 53 (3) (2015) 830–837.
201
[9] Piepenburg Olaf, et al., DNA detection using recombination proteins, PLoS Biol. 4 (7) (2006). [10] Lutz Sascha, et al., Microfluidic lab-on-a-foil for nucleic acid analysis based on isothermal recombinase polymerase amplification (RPA), Lab Chip 10 (7) (2010) 887–893. [11] H.T. Beernink, W.S. Morrical, The uvsY recombination protein of bacteriophage T4 forms hexamers in the presence and absence of single-stranded DNA, Biochemistry 37 (16) (2016) 5673–5681. [12] Shamoo Yousif, et al., Crystal structure of a replication fork single-stranded DNA binding protein (T4 gp32) complexed to DNA, Nature 376 (6538) (1995) 362–366. [13] Piepenburg Olaf, Niall A. Armes Recombinase polymerase amplification reagents and kits. U.S. Patent Application No. 13/375, 264. [14] S.C. Kowalczykowski, et al., On the thermodynamics and kinetics of the cooperative binding of bacteriophage T4-coded gene 32 (helix destabilizing) protein to nucleic acid lattices, Biophys. J. 32 (1) (1980) 403. [15] Jie Liu, Christopher L. Berger, Scott W. Morrical, Kinetics of presynaptic filament assembly in the presence of single-stranded DNA binding protein and recombination mediator protein, Biochemistry 52 (45) (2013) 7878–7889. [16] J.N. Farb, S.W. Morrical, Role of allosteric switch residue histidine 195 in maintaining active-site asymmetry in presynaptic filaments of bacteriophage T4 UvsX recombinase, J. Mol. Biol. 385 (2) (2009) 393–404. [17] R.D. Kuchta, et al., Kinetic mechanism of DNA polymerase I (Klenow), Biochemistry 26 (25) (1987) 8410–8417. [18] Mekbib Astatke, Nigel D. Grindley, Catherine M. Joyce, How E. coli DNA polymerase I (Klenow fragment) distinguishes between deoxy-and dideoxynucleotides, J. Mol. Biol. 278 (1) (1998) 147–165. [19] Xiao Ming, et al., Role of excess inorganic pyrophosphate in primer-extension genotyping assays, Genome Res. 14 (9) (2004) 1749–1755.
Contents lists available at ScienceDirect
Biochemical Engineering Journal journal homepage: www.elsevier.com/locate/bej
FA mathematical model of recombinase polymerase amplification under continuously stirred conditions Clint Moody a , Heather Newell a , Hendrik Viljoen a,b,∗ a b
Department of Chemical and Biomolecular Engineering, University of Nebraska-Lincoln, Lincoln, NE 68588, USA Preclinical Drug Development Platform, North-west University, Potchefstroom 2520, South Africa
a r t i c l e
i n f o
Article history: Received 18 January 2016 Received in revised form 18 March 2016 Accepted 19 April 2016 Available online 23 April 2016 Keywords: Recombinase Polymerase Amplification (RPA) Isothermal PCR Couette flow Mathematical model Disease diagnostics
a b s t r a c t Growing interest surrounds isothermal PCR techniques which have great potential for miniaturization for mobile diagnostics. Particularly promising, Recombinase Polymerase Amplification (RPA), combines this advantage of isothermal PCR with simplicity and rapid amplification. A mathematical model is presented of Recombinase Polymerase Amplification (RPA) and compared to experimental data. This model identifies the rate limiting steps in the chemical process, the effects of stirring, and insights in to using RPA for quantitative measurement of initial DNA concentration. Experiments are shown in which DNA amplification occurs under conditions of Couette flow and conditions of rotational turbulent flow. Hand mixing has been shown to dramatically shorten amplification times but introduces unpredictable variability. In some cases, this variability manifests itself as human error induced false negatives, a serious problem for all potential applications. Mechanical stirring demonstrates similarly short delay times while retaining high repeatability and reduces the potential for human error. Published by Elsevier B.V.
1. Introduction Polymerase Chain Reaction (PCR) has revolutionized molecular biology and has proven itself as one of the fastest and most specific methods of disease detection. In many places around the world, the sparsity of hospitals and clinics make developing point of care procedures a necessity. For PCR, the requirement for bulky and expensive thermocyclers remains one of the main obstacles for point of care use. To address this issue, many are turning towards isothermal PCR for the answer. A technique that does not require thermocycling to denature double stranded DNA but instead only requires a thermal bath lends itself much more readily to point of care diagnostics. Several isothermal techniques have emerged since the development of PCR. The most widely used are rolling-circle amplification, helicase-dependent amplification (HDA), Loop Mediated Isothermal Amplification (LAMP), and Recombinase Polymerase Amplification (RPA). Amplification times for Rolling circle generally take 65 min [1]. In addition to this 65 min process, it also requires significant preparation steps beforehand that can significantly drive up the total time
∗ Corresponding author. E-mail addresses: [email protected] (C. Moody), [email protected] (H. Viljoen). http://dx.doi.org/10.1016/j.bej.2016.04.017 1369-703X/Published by Elsevier B.V.
and overall complexity. The simplest and most efficient method for generating a suitable circular template, using a padlock probe, requires carefully designed oligomers and limits the amplification to small target regions [2]. Helicase-dependent amplification employs a helicase to unzip the double stranded DNA template to accommodate primer annealing. Helicase activity, however, requires bubbles to form naturally at AT rich regions [3]. The high GC content in Mycobacterium and other bacteria, in many cases, would inhibit this technique. Another drawback of HDA is that optimal amplification only occurs for short sequences of 80–120 bp [4]. Loop Mediated Isothermal Amplification (LAMP) has the benefit of being a simple method to perform from a kit but is deceptively complicated. It requires 4 primers that need to bind at 6 locations. The DNA regions around the primers must form loops without hairpins or other secondary structure so primer design is no trivial task [5]. Even though it does not require a high rate of thermocycling, it does require an initial denaturation step at 95◦ Celsius before amplifying at 65◦ C [6]. Amplifications with this method generally require 30–60 min [6]. A newly developed method called Recombinase Polymerase Amplification (RPA) shows great promise especially for pointof-care diagnostics. RPA uses recombinase to insert primers in to double stranded DNA rather than denaturing DNA using high temperature cycling [7]. This allows RPA to operate at a con-
194
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Nomenclature G m R FG n H3 PO4 FR FR FnR FnR D FD P PFD B np Sf t1 2
k1 k1f
Single stranded DNA binding protein (Gp32) Number of Gp32 binding sites on a primer Recombinase complex (UvsX.6*UvsY) Forward primer/Gp32 complex Number of recombinase binding sites on a primer Inorganic phosphate Unstable forward primer/recombinase complex Stable forward primer/recombinase complex Forward primer complexed with unstable filament of n recombinase molecules Forward primer complexed with unstable filament of n recombinase molecules DNA template Forward primer/DNA complex Polymerase Polymerase/forward primer/DNA complex Number of base pairs in template Number of base pairs in primer Empirical scaling factor to modify diffusion limited rate constants Time required for DNA concentration to reach 50% of saturation Equilibrium constant for Eq. (1) 105 M −1 Forward rate constant for Eq.
108 1/ (M × s)
5×
(1)
k3f
Reverse rate constant for Eq. (1) 5 × 103 1/s Equilibrium constant for Eq. (2a) 68 × M −1 Forward rate constant for Eq. (2a) 108 1/ (M × s) Reverse rate constant for Eq. (2a) 1.471/s Rate constant for Eq. (2b) 47 × 10−3 1/s Equilibrium constant for Eq. (2c) 3 × 106 M −1 Forward rate constant for Eq. (2c) 108 1/ (M × s) Reverse rate constant for Eq. (2c) 33 1/s Rate constant for Eq. (2d) 4.6 × 10−3 1/s The Michaelis constant for Reactions (3) and (4a) 20.35 × 10−6 M Forward rate constant for Eq. (3) 108 1/ (M × s)
k3r
Reverse rate constant for Eq. (3) 59.37 1/s
k1r keq2a keq2af keq2ar k2b keq2c keq2cf keq2cr k2d KM3a
k4acat k4bcat k5f
Forward rate constant for Eq. (4a) 4.22 1/s Forward rate constant for Eq. (4b) 8.32 1/s Forward rate constant for Eq. (5)
107 1/ (M × s)
2. Materials and methods 2.1. Amplicon design for real-time results 1.2 ×
k5r
Reverse rate constant for Eq. (5) 0.06 1/s
k6f
Forward rate constant for Eq. (6) 87 1/s
k7
Forward rate constant for Eq. (7) 4.1 1/ (M × s)
k8
Forward rate constant for Eq. (8) 1.13 1/ (M × s)
time PCR [7]. While many detailed models have been presented for the other isothermal PCR methods mentioned, RPA differs from them in some key ways that warrant special attention. To compare the model with data, the experimental methodology must minimize the effect of poorly repeatable inputs. One example of this would be hand mixing. The current methodology recommended by the supplier of RPA technology, Twistamp, requires two hand mixing steps. Our experience has shown that these mixing steps introduce a significant amount of variability in results. The tubes must be mixed both at the start of the process and after 4 min. The removal of the tube for mixing at 4 min creates a gap in data collection at a crucial point in the amplification process partly defeating the purpose of real time detection. To make matters worse, if the tubes are insufficiently mixed at 4 min, the samples will not amplify or amplification will be stunted. Another issue arises when, during the mixing step, droplets stick to the lids of the tubes artificially reducing the signal levels. These sources of human error make having a skilled operator imperative and confounds efforts to compare results obtained by different investigators. While these mixing steps are error prone and introduce significant randomness, without them the process takes several times as long to give a result. In this paper, we will demonstrate an automated stirring process that eliminates these sources of human error and simplifies the process to make it more easily employed in field applications. Automated stirring has been demonstrated before using lab-on-a-chip technology [10], but our method is more accessible to a wide range of users without buying specialty disposable instrumentation. Initially, our procedure used an ordinary benchtop drill press to stir the reaction mixture. The equipment was then scaled down to use instead a small RC motor. The only part of the set-up that needed replacement after each test was the drill bit/stir rod. Not only does this motor take up negligible space (about as much as a 1.5 ml tube) but its low cost also makes an array of 8, one for each tube, quite feasible even in low budget applications. This stirring method, by removing the arbitrary element of hand shaking and moving the reaction in to a kinetically controlled regime, made the reaction more conducive to mathematical modeling. We herein present, to our knowledge, the first mathematical model specifically addressing the unique considerations of RPA.
stant temperature of around degrees celsius—much lower than other isothermal techniques like LAMP [8]. Using the tubescanner from TwistDx, the instrumentation costs for RPA are significantly cheaper than those for either PCR or LAMP and include fluorometers for real time detection [8]. It only requires two primers that, with a length of 30–35 base pairs, are only modestly larger than PCR primers [9]. RPA also boasts shorter amplification times than LAMP. In combination with real time fluorescent measurements, a result can be determined in approximately 10 min [7]. Remarkably, the rapid rate of amplification does not come at the expense of reduced sensitivity. RPA can amplify as few as 10 initial copies of template—a detection threshold similar to that of traditional real
DNA was extracted from Mycobacterium smegmatis and Escherichia coli (as a negative control) using the phenol-chloroform method. The amplicon for the detection of the Mycobacterium Smegmatis DNA was designed using the sequence IS1096. The RPA amplicon was 116 base pairs long. SYBR green was used instead of designing a specific probe. The forward and reverse primers were: IS1096F: 5 -CTCATCGAACATTCCCGCGAACACGTTCCGACCAG-3 IS1096R: 5 -CTTGACGGTGTAGAGACGATCAGCTGCTTTCGC-3 2.2. RPA conditions and detection RPA was performed using the TwistAmp Basic kit, having a 50 l volume. Reaction mixtures were formulated according to TwistDx Ltd. recommendations with one of notable change. Instead of using probes for detection, 5× concentrated SYBR Green was added to each reaction mixture. The DNA concentration was standardized at 0.1 ng/l for all tests. Once all components were added, the tubes were vortexed and centrifuged; a solution of 14 mM Mg acetate was then pipetted into each reaction tube which was then immediately put into the tubescanner device. The reaction components were
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
supplied by TwistDx Ltd., Cambridge, United Kingdom. The fluorescent signal is in the FAM channel, 470 nm excitation and 520 detection. Each experiment was run for 20 min at 39◦ C. The realtime readout was courtesy of Twista Studio software, version 2.6.0 (TwistDx, Ltd.) which interpreted fluorescent data into graphs of Intensity (mV) vs. time (s). 2.3. Experimental set-up For all reactions stirred at or below 3050RPM, RPA reaction mixtures were stirred with a PerformaX 10” benchtop drill press using 1/16th in. stainless steel drill bits. Following each experiment, the drill bit was replaced to avoid contamination. Experiments using 1/16th in. diameter polypropylene rods (Data not shown) gave similar results but required higher rotation speeds. In order to allow the drill access to the reaction mixture, the lid of the TwistA system had to be removed. Since this exposed the optical detection system to ambient lighting, all experiments were performed with the lights out. For our highest stirring speed test, at 7050RPM, we used an RC motor. The exact model of motor was the Ares AZSH1313 Ares Ethos QX130 Motor w/Pinion Gear and Wire Leads, Clockwise Rotation. A brass coupler and two set screws were used to mount the drill bits to the RC motor. This particular motor becomes unreliable at lower speeds and thus could not be used in the tests below 7050RPM. In many applications, high speed stirring provides the best results. Only when used with oligomer probes do we see improved performance with slower stirring speeds. For these applications, other similar motors with slower optimal speeds could be substituted. 3. Theory PCR employs thermal cycling to denature double stranded DNA and facilitate primer annealing. Recombinase polymerase amplification, instead, uses a recombinase complex to denature double stranded DNA so that primers can bind. This recombinase complex consists of an UvsY hexamer bound to UvsX. In our model, we assume that UvsY forms a stable hexamer and binds to UvsX prior to the start of the test. Under normal RPA conditions, this is a safe assumption because the protein pellet contains both proteins and the UvsY hexamer is stable with a small dissociation constant [11]. Another essential element of the RPA pellet, a single stranded DNA binding protein called Gp32, helps remove secondary structure from primers to facilitate the binding of recombinase [12]. We assume that Gp32 binding must proceed recombinase binding. Thus the first reaction treated by our model involves primers binding reversibly with Gp32, represented simply by G. m × G + F FG
(1)
For the sake of simplicity, we assume perfect symmetry between the reactions involving forward and reverse primers. This allows us to reduce the number of differential equations in the model by simply multiply the appropriate terms by a factor of two. The RPA pellet composition as a patented technology is publically available and, assuming a 50 ul reaction, we calculate a 9.4 × 10−7 Mconcentration of Gp32 [13]. The published value for the equilibrium constant of this reaction, k1 , is 105 M −1 [14]. The dissociation rate k1r is5 × 103 1/s [14]. The binding site for Gp32 is 7 base pairs long so for a primer 32 base pairs long m is roughly 4 [14]. Next, the recombinase complex, R, displaces Gp32 from the primer in a two-step process [15]. First, the process requires a nucleation step. R + FG FGR
FGR → FGR + G
(2a) (2b)
195
The complex binds to the primer to become the unstable complex FGR but only becomes stable when it displaces Gp32 to form FGR [15]. While UvsX uses ATP hydrolysis to drive this process, we assume that ATP levels are constant and not limiting. The presence of phosphocreatine kinase and phosphocreatine in the pellet acting to regenerate ATP makes this assumption justified. Once one UvsX.UvsY complex has bound and achieved stability, nucleation is complete and the process can proceed to the growth phase. FGR + (N − 1) × FGnR
(2c)
FGnR → FnR(m − 1) × G
(2d)
Each UvsX monomer binds to 4 base pairs of DNA so the number of monomers consumed, n, is approximately 8 for primers of length, np = 32 [16]. Liu et al. lists the kinetic parameters and equilibrium constants for these four reactions above [15]. This paper presents alternative values using either hyperbolic or global methods. The global parameters better fit the data so they were used in all four reactions. The value of the equilibrium constant of Reaction (2a), keq2a ,as determined by global methods is 68 × M −1 [15]. The value of the forward rate constant of reaction (2b), k2b , is 47 × 10−3 1/s [15]. The value of keq2c , another equilibrium constant this time for Reaction (2c), for kinetically controlled reactions is 3 × 106 M −1 [15]. The rate constant for the complete dissociation of Gp32, k2d , is 4.6 × 10−3 1/s [15]. Based on the pellet composition the concentration of UvsX is 5.9 × 10−6 M [13]. Once recombinase has bound a primer, it seeks out homology in the DNA template. It then inserts the primer in to the DNA and dissociates from both by hydrolyzing ATP. We originally modeled this reaction using Michaelis-Menten kinetics as in Farb et al. We quickly discovered, however, that the underlying assumptions of the Michaelis-Menten kinetics poorly described this system. So we instead used the more fundamental two reaction formulation. FnR + D + 2 × m × G FnRD
(3)
FnRD + n × ATP → FD + n × R + n × AMP + n × PPi
(4a)
In Reaction (3), the recombinase filament escorts the primer to the region of homology in the DNA template, here represented by D. Next, hydrolysis of ATP drives the dissociation of the recombinase from the primer/DNA complex, FD. Michaelis-Menten kinetics assumes that the equilibrium reaction has reached steady state. This assumption in particular does not apply during the transient when diffusion rates are slow and substrate concentrations are low. We used the constants k4acat = 4.22 1/s and, for the unstirred reaction, KM3a = 20.35 × 10−6 M [16]. Along with the assumption that the forward reaction rate k3f = 108 1/ (M × s), corresponding to the diffusion limited case, we were able to derive the rate constants from the Michaelis-Menten parameters. This gives us a reverse reaction rate, k3r, of 59.37 (1/s). The binding of GP32 to the opposing DNA strand helps stabilize the complex formed in Reaction (3). Unlike RecA and other members of the recA family, UvsX has two modes of ATP hydrolysis. It can hydrolyze ATP to either AMP or ADP [16]. As side products, Reaction (4a) produces pyrophosphate (PPi) and Reaction (4b) produces inorganic phosphate (H3 PO4 ).To account for this secondary mode of ATP hydrolysis we added the following reaction. FnRD + n × ATP → FD + n × ADP + H 3 PO4
(4b)
We used the forward rate constant k4bcat = 8.32 1/s for this reaction [16]. RPA employs strand displacing Sau Polymerase (a polymerase derived from Staphylococcus aureus denoted here as P) to extend primers [7]. Polymerase can bind reversibly to DNA primer complexes. We assume binding rates comparable to those
196
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 1. Overall reaction scheme diagram of RPA.
of DNA polymerase I which are well established as k5f = 1.2 ×
107 1/ (M × s) and k5r = 0.06 1/s [17]. FD + P PFD
(5) −6
The concentration of Polymerase is 1.34 × 10 M [13]. Once the tertiary complex has been formed polymerase can extend the primer to form a new strand of DNA using nucleotides (dNTP’s). The factors of two in this equation indicate that a simultaneous extension of the reverse primer generates a second strand of DNA. Again we must assume a rate constant of similar magnitude as E. coli Polymerase I which has been measured to be approximately k6f = 87 1/s [18]. PFD + 2 × B × dNTP → P + 2 × B × PPi + 2 × DNA
(6)
If pyrophosphate exists in sufficiently high quantities, it can drive the reverse reaction [19]. The reverse reaction rate was calculated from the forward reaction rate and the Gibbs free energy change of the extension process. Using the Gibbs free energy reported in Kuchta et al., the calculated the equilibrium constant is approximately 77 [17]. Thus k8 = k6f /77. Since this degradation is proportional to the length of DNA k7 = k8 × B/32. This process had little effect on the model outcome under the experimental conditions examined but may be significant for some reactions. FD + np × PPi → D + np × dNTP
(7)
D + B × PPi → B × dNTP + FD
(8)
This paper focuses mainly on the TwistDx Basic and FPG kits. These pellet formulations do not include exonuclease. The equations, however, can be easily modified to accommodate its presence. Reactions for exonuclease activity would take the same form as (7) and (8) (Fig. 1). Before presenting the rate equations for RPA, it is important to note that the pieces of the model came from multiple papers and we needed to make the different formulations compatible. We define DNA concentration on a per copy basis but our sources define it based on the number of binding sites or the number of base pairs. The rate constants for the Gp32 equilibrium used Gp32 binding sites as their reference. The rate and equilibrium constants for the recombinase reactions were calculated per recombinase binding site. Thus we need correction factors for our per primer calculations. The rate constants for primer insertion in to dsDNA were defined based on a per base pair concentration of DNA rather than the copy number so correction factors account for this difference in methodology (Fig. 2). The rate equations we use are: d [R] = −2 × keq2af × [R] × (n2 × [FG]) + 2 × keq2ar dt
× FGR − 2 × keq2cf × [R] × FGR + 2 × keq2cr × FGnR
+ 2 × k4acat × [FnRD] + 2 × k4bcat ∗ [FnRD]
Fig. 2. Model predictions for DNA concentration compared with three replicates of experimental data with a stirring speed of 7050RPM.
d ([FG]) = 1/m × k1f × [G] × [F] − 1/m × k1r × [FG] − keq2af dt
× [R] × (n × [FG]) + 1/n × keq2ar × FGR )
d [FGR ] = keq2af × [R] × (n × [FG]) − 1/n × keq2ar dt
× FGR − k2b × FGR
d [FGR] 1 = k2b × FGR − × keq2cf × [R] × [FGR] n−1 dt +
1 × keq2cr × FGnR n−1
1 d [FGnR ] × keq2cf × [R] × [FGR] = n−1 dt −
1 n−1
× keq2cr × FGnR − k2d × FGnR
d [FnR] = k2d × FGnR + dt ×k3f × [FnR] × [D]
1 n
× k3r × [FnRD] −
1 n
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 3. Model predictions for the concentration of primers forming an unstable complex with one molecule of recombinase. d[FnRD] dt
=
1 n
k4acat × [FnRD] − d[D] dt
=−
1 n
× k3f × [FnR] × [D] −
1 n
1 n
× k3r × [FnRD] −
n
Fig. 4. Model predictions for the concentration of primers forming a stable complex with one molecule of recombinase.
×
× k4bcat × [FnRD]
× k3f × [FnR] × [D] +
1 n
× k3r × [FnRD] + 2 ×
k6f1 × [PFD] + 2 × k7 × [PPi] 1 × [FD] − 2 × k8 × [PPi] ∗ [D] n
1
197
× k4acat × [FnRD] +
n
d[FD] dt
=
× k4bcat × [FnRD] − k5f × [P] ×
[FD] + k5r × [PFD] − k7 × [PPi] × [FD] + k8 × [PPi] × [D] −k5f × [P] × [FD] + k5r × [PFD] + k6f × [PFD]
d[P] dt
=
d [PFD] = k5f × [P] × [FD] − k5r × [PFD] − k6f × [PFD] dt d [dNTPs] = −2 × B × k6f × [PFD] + 2 × B × k8 × [PPi] dt × [D] + 2 × np × k7 × [PPi] × [FD]
d [PPi] = 2 × k4acat × [FnRD] + 2 × B × k6f × [PFD] − 2 × B dt ×k8 × [PPi] × [D] − 2 × np × k7 × [PPi] × [FD]
4. Results and discussion Using the model described above, we compared with experimental results and found reasonable agreement with data taken at 7050RPM. This best corresponds to one of the fundamental assumptions of this model: namely that the fluid is well mixed and all reactions are kinetically controlled. The Reynold’s number of this rotational flow is roughly 5000 so turbulence will dramatically enhance the mixing of the reactants. Without consistent mixing to minimize diffusion effects, there would be no foundation for this model. Later, we will relax the kinetically controlled assumption to describe the effect of slower stirring rates but first we must characterize the simplest form of the model. To understand the time evolution of the RPA reaction, we simulated the differential equations using a second order Euler method. Our simulations showed that the primers and Gp32 bind almost completely and instantaneously. So, for appropriate Gp32 concentrations, one could safely assume that all primers exist as this complex at t = 0s with negligible error. Likewise, FGR’ appears almost instantaneously (Fig. 3).
Fig. 5. Model predictions for the concentration of primers forming an unstable complex with a filament of n molecules of recombinase. In this case n = 8.
The first step to take any significant time (2b) produces peak levels of FGR at around 50s (Fig. 4). The production of FGnR’ (2c) generally tracks the concentration of FGR but introduces significant smoothing (Fig. 5). Fig. 6 shows that the rate limiting reaction for RPA is Reaction (2d) which suggests that the process could be improved by reducing the binding affinity of the single strand binding protein. The time evolution of the next three species produced (FnRD, FD, and PFD) appear qualitatively similar to the concentration of FnRD in Fig. 7. Under the conditions of our experiments, not only are they qualitatively similar, they roughly proportional. Estimating PFD as 0.475*FnRD, only results in a 4% error in FnRD across the peak in concentration (Fig. 8). Thus we could simplify the model by eliminating two differential equations while generating nearly identical results. This estimate is valid for a wide range of initial DNA concentrations so long as other variables such as polymerase concentration remain constant. For polymerase concentrations below optimized levels, FnRD and FD remain roughly proportional but limited polymerase concentrations control the production of PFD.
198
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 6. Model predictions for the concentration of primers forming a stable complex with a filament of n molecules of recombinase. In this case n = 8.
Fig. 8. Model predictions showing that the concentration of polymerase/primer/DNA complex is roughly proportional to the concentration of primers being delivered to the DNA template.
Fig. 7. Model predictions for the concentration of recombinase delivering primers to the DNA template.
The key reaction in the process is (3). The limiting reactant in this process controls the rate of DNA production. RPA, for the conditions investigated, occurs in three distinct stages. During the first stage, the assembly of the required FnR occurs. Once these primer/recombinase complexes form and far exceed the number of available DNA molecules, the second stage begins. The second stage shares much in common with the various methods of isothermal PCR. During this stage, the limited supply of DNA determines the rate and DNA concentration follows a roughly exponential curve. Once DNA reaches high enough levels, FnR once again becomes rate limiting. Without measurements of FnR, it is impossible to see the transition between the first and second stages based on a DNA plot alone (except in the unusual circumstance that DNA concentrations begin at detectable levels). The transition between the second and third stages, however, can be seen as a rounded corner in the graph which ends the exponential growth phase. These three stages result in behaviour that differs qualitatively from traditional PCR. Since FnR production occurs independently
Fig. 9. Model predictions of the effect of varying polymerase concentrations on the rate of DNA production.
from reactions after (3), any changes after that point will only affect phase 2. If some quantity such as polymerase concentration changes, the graphs continue to converge at the beginning of stage 3. Changes made that effect reactions before (3), however, affect all three stages. The major contribution to the delay time, under our experimental conditions, is the extension step of the recombinase filament assembly. This contrasts sharply with the time evolution in other forms of isothermal PCR. While the graphs may look qualitatively similar, the delay time in RPA is much less sensitive to the initial DNA concentration. As can be seen from Fig. 9, the polymerase concentration has been optimized for the commercial pellets and any further increases in concentration or improvements in the polymerase processivity would have almost no effect. Decreasing polymerase concentration, on the other hand, can significantly increase the
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Fig. 10. Model predictions of the effect of different initial DNA concentrations on DNA production.
199
Fig. 11. Model predictions showing increased sensitivity to DNA concentration when Polymerase concentration is reduced by a factor of 100.
delay time. This second effect may find unexpected usefulness in quantitative studies of initial DNA concentration. As has been shown previously, using quantitative PCR techniques one can arrive at a reasonable estimate of the initial concentration of DNA simply extrapolating backwards assuming exponential growth. This same method naively applied to RPA will result in nonsensical estimations orders of magnitude smaller that the true initial value. For the current optimized pellet we must instead use the shifted equation: t 1 = −14.5 × log10 (D) + 16.5 2
where t 1/2 is the time required for the signal to reach ½ its maximum value. The time shift accounts for the time required to assemble the recombinase filament on primers. Fig. 10 shows the effect of different initial DNA concentrations and one immediately notices that in its current form RPA is not terribly conducive to quantitative studies because such small differences in delay time can easily fall below the experimental variability. Reducing polymerase concentration or employing a less processive polymerase could potentially resolve this problem. Both methods would slow DNA extension while leaving the RPA specific reactions relatively unaffected. Since commercial pellets contain both the polymerase and the recombinase, it may be advisable to reformulate the pellet composition for quantitative RPA. Simply reducing the quantity of the pellet would reduce the recombinase levels below the concentration needed for the reaction to proceed. To demonstrate this spreading effect DNA dilutions were simulated at the current Polymerase concentration as well as at a hundredth of that concentration (Fig. 11). For the optimized polymerase concentration the delay time to half maximum is given by the formula above. For the reduced polymerase case, the formula is instead:
Fig. 12. Model predictions for the effect of varying recombinase concentrations on the rate of DNA amplification.
dramatic effect on the overall amplification. Amplification drops to nearly zero at 1/3 of the required recombinase. At these low concentrations, primers must compete for scarce recombinase complexes and complete filament assembly becomes exceedingly rare. Thus the entire process becomes stalled between the nucleation step (2b) and the growth step (2c). 4.1. Experimental stirring
t 1 = −33.2 × log10 (D) + 98.4 2
The spreading increased by a factor of 2.3 while the delay has only increased by a factor of 1.7. Using the right pellet formulation, this method could potentially generate quantitative results of any desired resolution. The effect of varying concentrations of recombinase was also simulated. From Fig. 12, we can see the insensitivity of DNA amplification to excess recombinase. Once recombinase levels become insufficient, however, small differences in concentration have a
Initially, a concern with mechanical stirring of the reaction mix was that the drill bit would interfere with the optical signal. Fortunately, experiments showed that a 1/16th in. drill bit submerged 3.5 mm in to the fluid had only a small effect on the signal produced by a 50 l reaction volume. For larger diameter drill bits, a constant offset was observed which could easily be subtracted from the signal to reproduce the expected results. While TwistDx recommends using specialty DNA oligomer probes, we did not observe any issues with using SYBR green either in terms of amplification
200
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
time for all stirring speeds investigated. The equation of the best fit line with an R2 value of 0.51 was: t 1 = −0.0082 × RPM + 273 2
Fig. 13. Experimental results for different mixing rates averaged over 3 replicates.
Tests done with stirring speeds below 570RPM amplified but their time evolution displayed poor repeatability (Data not shown). This poor repeatability was especially pronounced in experiments conducted with no mixing at all. While the rate of reaction can be dramatically improved by hand mixing, this introduces significant problems for repeatability as well. Aside from the inherently arbitrary nature of hand mixing, care must also be taken so that droplets do not stick to the lid of the tubes. These droplets contain DNA that will be hidden from the detector and thus can greatly diminish the observed fluorescence. Most significantly, we observed a higher rate of false negatives with hand mixing than with drill stirring. For any application, false negatives should be kept to a minimum and any attempt at quantifying the initial DNA concentration using RPA will require a mixing method with highly repeatable outcomes. Now we are ready to discuss modifications to the model that incorporate different stirring speeds. The forward reaction rates of (2a), (2c), and (3) are all diffusion limited. Thus we expect the corresponding reaction rate constants to decrease as stirring speed decreases. To model this we introduce a scaling factor, Sf, to reduce the forward rate constant for each of these reactions. Based on our experimental data, the appropriate Sf can be roughly calculated as:
Sf = 8.05 × 10−5 × RPM + 0.347 To use this scaling factor, simply replace keq2af, keq2cf, and k3f with Sf*keq2af, Sf*keqcf, and Sf*k3f respectively. We determined the scaling factor for each experimental condition by generating a simulation with the same delay time as our data collected at the different stirring speeds. This represents merely a simple empirical formulation which may prove useful to the experimentalist. The underlying chemical phenomenon cannot be determined from our limited data but clearly involves the enhancement of diffusion controlled Reactions—(2a), (2c), and (3) being the most significant. Since this data was taken at Reynolds numbers ranging from roughly 400–2200 (all within the regime of laminar flow), it should not be extrapolated above 3200RPM at which point the flow begins to transition to turbulent flow. Fig. 14. Experimental results showing the relationship between the time required for the reaction to reach half its maximum value and the stirring speed in RPM.
5. Conclusions inhibition or increased noise. Since probes must be tailored specifically for each application which represents a significant cost, and must be optimized empirically, we recommend using SYBR green in applications that do not require multiplexing and exhibit minimal nonspecific amplification. Using E. coli DNA, we showed that no detectable non-specific amplification occurred in one hour of incubation. So, unless otherwise noted, all experiments were conducted with SYBR Green however similar results were obtained with the recommended oligomer probes. The only significant difference we observed between the probes and SYBR Green was that when drill speeds exceeded 900RPM it compromised the probes’ ability to fluoresce. Probes appeared to work normally at drill speeds of 570RPM so the system will work with probes given the proper choice of drill speed. Increasing rotation speed increases the rate of amplification (Fig. 13). The relationship between stirring speed measured in RPM and the delay time of the amplification is roughly linear (Fig. 14). We continue to see a small but consistent improvement in delay
1. The equations listed represent to our knowledge the first mathematical model presented for Recombinase Polymerase Amplification 2. Consistent quantitative experimental results can be achieved using mechanical stirring eliminating human variability as well as making the process more reliable. 3. Stirring significantly increases the rate of amplification. 4. Insufficient recombinase dramatically hinders or even stalls the reaction. 5. Amplification can be understood as occurring in 3 stages: the first and third limited by FnR while only the second is limited by DNA concentration. The presence of the first stage makes estimating the initial DNA concentration different than the methods used in traditional quantitative PCR. Exponential growth models can underestimate initial DNA concentration by several orders of magnitude. 6. The model gives insights of value to anyone attempting to use RPA to quantify DNA concentration.
C. Moody et al. / Biochemical Engineering Journal 112 (2016) 193–201
Acknowledgements We thank D.K. Zinniel for providing DNA samples from Mycobacterium Smegmatis and E. coli. References [1] Chia-Chen Chang, et al., Diagnostic devices for isothermal nucleic acid amplification, Sensors 12 (6) (2012) 8319–8337. [2] Jia Li, Joanne Macdonald, Advances in isothermal amplification: novel strategies inspired by biological processes, Biosens. Bioelectron. 64 (2015) 196–211. [3] Schaerli Yolanda, et al., Isothermal DNA amplification using the T4 replisome: circular nicking endonuclease-dependent amplification and primase-based whole-genome amplification, Nucleic Acids Res. 38 (22) (2010) e201. [4] Kong Huimin, Myriam Vincent, Yan Xu, Helicase dependent amplification of nucleic acids. U.S. Patent No. 7, 282, 328. 16 Oct. (2007). [5] Notomi Tsugunori, et al., Loop-mediated isothermal amplification of DNA, Nucleic Acids Res. 28 (12) (2000) e63. [6] Tomita Norihiro, et al., Loop-mediated isothermal amplification (LAMP) of gene sequences and simple visual detection of products, Nat. Protoc. 3 (5) (2008) 877–882. [7] Euler Milena, et al., Recombinase polymerase amplification assay for rapid detection of Rift Valley fever virus, J. Clin. Virol. 54 (4) (2012) 308–312. [8] Teoh Boon-Teong, et al., Early detection of dengue virus by use of reverse transcription-recombinase polymerase amplification, J. Clin. Microbiol. 53 (3) (2015) 830–837.
201
[9] Piepenburg Olaf, et al., DNA detection using recombination proteins, PLoS Biol. 4 (7) (2006). [10] Lutz Sascha, et al., Microfluidic lab-on-a-foil for nucleic acid analysis based on isothermal recombinase polymerase amplification (RPA), Lab Chip 10 (7) (2010) 887–893. [11] H.T. Beernink, W.S. Morrical, The uvsY recombination protein of bacteriophage T4 forms hexamers in the presence and absence of single-stranded DNA, Biochemistry 37 (16) (2016) 5673–5681. [12] Shamoo Yousif, et al., Crystal structure of a replication fork single-stranded DNA binding protein (T4 gp32) complexed to DNA, Nature 376 (6538) (1995) 362–366. [13] Piepenburg Olaf, Niall A. Armes Recombinase polymerase amplification reagents and kits. U.S. Patent Application No. 13/375, 264. [14] S.C. Kowalczykowski, et al., On the thermodynamics and kinetics of the cooperative binding of bacteriophage T4-coded gene 32 (helix destabilizing) protein to nucleic acid lattices, Biophys. J. 32 (1) (1980) 403. [15] Jie Liu, Christopher L. Berger, Scott W. Morrical, Kinetics of presynaptic filament assembly in the presence of single-stranded DNA binding protein and recombination mediator protein, Biochemistry 52 (45) (2013) 7878–7889. [16] J.N. Farb, S.W. Morrical, Role of allosteric switch residue histidine 195 in maintaining active-site asymmetry in presynaptic filaments of bacteriophage T4 UvsX recombinase, J. Mol. Biol. 385 (2) (2009) 393–404. [17] R.D. Kuchta, et al., Kinetic mechanism of DNA polymerase I (Klenow), Biochemistry 26 (25) (1987) 8410–8417. [18] Mekbib Astatke, Nigel D. Grindley, Catherine M. Joyce, How E. coli DNA polymerase I (Klenow fragment) distinguishes between deoxy-and dideoxynucleotides, J. Mol. Biol. 278 (1) (1998) 147–165. [19] Xiao Ming, et al., Role of excess inorganic pyrophosphate in primer-extension genotyping assays, Genome Res. 14 (9) (2004) 1749–1755.