Reaction Kinetics Simulation

Introduction

Reaction kinetics simulation is a mathematical modeling method aimed at predicting the trends of concentration of reactants or products with certain parameters by constructing a system with ordinary differential equations. We simulated the cleavage process of CRISPR/Cas13a system with the Simbiology toolbox in MATLAB^[1]. By applying the Michaelis-Menten equation, we simulated the cleavage of reporter RNA by our engineered CRISPR/Cas13a system under different conditions.

Simulation of Michaelis-Menten equation

Once miR-21-5p binds to crRNA, it activates the cleavage activity of the Cas13a/crRNA complex, which cleaves not only miR-21-5p but also other RNA. Therefore, the Cas13a/crRNA/miRNA ternary complex can be treated as an enzyme with RNA cleavage activity. We use a fluorophore (Carboxyfluorescein, FAM) and a quenching group (Black hole quencher, BHQ1) dual-labeled 6U oligonucleotides as a reporter RNA to evaluate the cleavage efficiency of the Cas13a/crRNA/miRNA ternary complex, since cleavage of the reporter RNA results in the release of fluorescence quenched by BHQ1, which can be measured easily. Such process can be simplified as the reaction below.

The positive reaction rate constant for the binding reaction between the Cas13a/crRNA/miRNA ternary complex and the reporter RNA is k1, the reverse reaction rate constant is k2, and the reaction rate constant for fluorescence generated by cleavage is k3. We can use Michaelis-Menten equation to study the rate of cleavage activities^[2].

\[ v_0 = \dfrac{dF}{dt} = \dfrac{V_m[S]}{K_m+[S]} \]

\(v_0\) represents reaction rate, which is experimentally expressed as the rate of increase in fluorescence intensity, and [S] represents the concentration of the substrate, which is the reporter RNA. In actual reactions, fluorescence gradually degrades. Therefore, the equation can be rewritten as below.

In this equation, the \(k_{cleavage}\) represents the rate of the enzyme activity, which is the cleavage activity of Cas13a/crRNA/miRNA complex in our case. \(k_{cleavage} = k3\), \( K_m=\frac{k2+k3}{k1}\). We set the initial concentration of the ternary complex as 20 nM and the initial concentration of the reporter RNA as 1000 nM. The parameter settings are as follows.

\[ k1 = 10 nanomole^{-1} {hour}^{-1} \]

\[ k2 = 0.1 {hour}^{-1} \]

\[ k3 = 1000 {hour}^{-1} \]

\[ r = 0.5{minute}^{-1} \]

By changing the initial concentration of substrate, enzyme, and the values of different parameters, the following simulation results can be obtained.

Figure 1: Curves of fluorescence production with different association rate constants

Figure 2: Curves of fluorescence production with different decomposition rates

Figure 3: Curves of fluorescence production with different degradation constants

Figure 4: Curves of fluorescence production with different reporter RNA concentrations

It can be seen from (1) that the larger the forward reaction rate constants (\(k_1, k_3\)), the faster the fluorescence reaches its peak. The larger the degradation constant, the faster the fluorescence reaches its peak. Considering the degradation of fluorescence, the fluorescence intensity ultimately tends to zero. Since the main difference in actual reactions is the miRNA concentration, and under physiological conditions miRNA will not be saturated. Therefore we use miRNA concentration to represent the concentration of the ternary complex in our model. Setting the r value to 0.1 and changing the initial concentration of the ternary complex, the following curves are obtained.

Figure 5: Curves of fluorescence production with different ternary complex concentrations

Considering fluorescence degradation, there is an approximate linear relationship between fluorescence intensity and miRNA concentration within a large range of miRNA concentration (Figure 5).

Experimental validation

An effective model should fit the results of experiments, so we validated the model through experimental data. Fluorescence was measured from the actual experiments (Figure 6 & 7). In figure 6, the concentrations of crRNA, LwaCas13a, and miRNA were 20 nM, 100 nM, 50 nM respectively. In figure 7, the concentrations of crRNA, LwaCas13a, and reporter RNA were 20nM, 100nM and 100nM respectively.

Figure 6: Fluorescence intensity measured at different times with different concentrations of reporter RNA

Figure 7: Fluorescence intensity measured with different concentrations of miR-21-5p

In Figure 7, from 0 to 100 nm of miRNA, the fluorescence intensity increased with the increased concentration of miRNA, and it reached a peak at 100 nm of miRNA, indicating that the formation of Cas13a/crRNA/miRNA ternary complex were likely saturated by 100 nm of miRNA. Subsequently, as the concentration of miRNA increased to above 100 nm, the excess miRNA competed against reporter RNA as the substrate of the ternary complex, and suppressed the cleavage on reporter RNA, which decreased the fluorescence intensity.

In the increasing part of the fluorescence curve, r can be ignored, and taking the reciprocal on both sides of the equation yields the following equation:

\[ \dfrac{1}{v_0} = \dfrac{K_m}{V_m}\dfrac{1}{[reporter]}+\dfrac{1}{V_m} \]

Value before peak and fit it to obtain the following curve:

Figure 8: The fitting curve of Lineweaver-Burk

It can be seen that after using Lineweaver-Burk plot, the experiment results comply with the Michaelis-Menten equation, which is consistent with our predictions

Conclusion

Both experiment results and reaction kinetics simulation demonstrate that our CRISPR/Cas13a system can detect miR-21-5p at a large range of concentration. This feature of our CRISPR/Cas13a system is very important, since the concentration of miR-21-5p might vary a lot in different person.

References

[1] MATLAB and Simbiology Toolbox Release 2020b, The MathWorks, Inc., Natick, Massachusetts, United States.

[2] Ramachandran A, Santiago JG. CRISPR Enzyme Kinetics for Molecular Diagnostics. Anal Chem. 2021 May 25;93(20):7456-7464. doi: 10.1021/acs.analchem.1c00525.