Furin and the adaptive mutation of SARS-COV2: a computational framework.

SARS-2 virus has reached its most harmful mutated form and has damaged the world's economy, integrity, health system and peace to a limit. An open problem is to address the release of antibodies after the infection and after getting the individuals vaccinated against the virus. The viral fusion process is linked with the furin enzyme and the adaptation is linked with the mutation, called D614G mutation. The cell-protein studies are extremely challenging. We have developed a mathematical model to address the process at the cell-protein level and the delay is linked with this biological process. Genetic algorithm is used to approximate the parametric values. The mathematical model proposed during this research consists of virus concentration, the infected cells count at different stages and the effect of interferon. To improve the understanding of this model of SARS-CoV2 infection process, the action of interferon (IFN) is quantified using a variable for the non-linear mathematical model, that is based on a degradation parameter γ . This parameter is responsible for the delay in the dynamics of this viral action. We emphasize that this delay responds to the evasion by SARS-CoV2 via antagonizing IFN production, inhibiting IFN signaling and improving viral IFN resistance. We have provided videos to explain the modeling scheme.

Introduction

The impact of the CoViD-19 disease is devastating and scientific research seek to understand the mechanisms of infection to create vaccine with least side-effects and most promising results against the mutated virus. Scientists have reported that Arg-Arg-Ala-Arg (RRAR) is a cleavage site for furin enzyme. The furin action on S protein permit a faster activation and greatest rate of infection.. The furin site Is on SARS-COV2 but not in SARS-CoV or MERS-CoV. It is responsible for the high infection rates and transmission rates of SARS-CoV2. The presence of the this cleavage site is experimentally proven Walls et al. (2020) and the Activation of S requires proteolytic cleavage at two distinct sites: in the unique multibasic site motif RRAR, located between the S1 and S2 subunits, and within the S2 subunit (“S2”) located immediately upstream of the hydrophobic fusion peptide that is responsible for triggering virus-cell membrane fusion. This event, although not exclusive to SARS-CoV2, is important because it is absent on the viral antigens of the same viral family Coutard et al. (2020); Yu et al. (2021). It is important to consider this characteristic because it implies a greater speed of action of SARS-CoV2 than, for example, SARS-CoV, so much so that some adaptive mutations such as D614G seem to carry out structural changes that more expose the cut site for furin Korber et al. (2020).

In this article, we are focusing on the speed of action of SARS-2 inspired by the experimental study of Papa et’al. Papa et al. (2021). The exclusive action of the cutting site of furin, that is present in SARS-2 is still an open problem.

The work of Buonvino and Melino (2020) identifies viral evolution from the RaTG13 genotype and shows how, for SARS-CoV2, the acquisition of the furin cleavage site implies greater instability of the S protein. This is a very important factor for understanding the dynamics of the infection. For the infectious process to begin, some key enzymes play important role. For example, Furin” contributes to split the S protein into two subunits: S1 and S2. Therefore, it facilitate the fusion between the viral membrane and that of the host cell. SARS-CoV2 presents a further modified cleavage site for furin, i.e. in the amino acids of this site, a proline is added that changes the sequence and allows a strong bending of the structure leading to the introduction of three glycans O-linked that line the site itself. Furthermore, the furin-promotes infection capacity as well as the adaptive mutation. SARS-CoV2 virus has acquired a cleavage site for furin between S1 and S2 that appears to promote pathogenicity. This fact, in addition to enhancing the viral pathogenic aspect, also seems to be responsible for the speed of infection, especially in connection with the presence of an adaptive mutation “D614G”. “D614G mutation” neither increases S protein affinity for ACE2 nor makes viral particle more resistant to neutralization and that TMPRSS2 and Furin of all species studied can cleave the SARS-CoV-2 S glycoprotein in a similar way, provided that they are well conserved proteases among many species Brooke and Prischi (2020). For further details, please see video S1 and S2.

By taking into account these properties of Furin, we have hypothesized that the high rate of infection occurs when there is the presence of this adaptive mutation together with a structural adaptation of the clevage site of the furin on the viral protein S.

It is highly desired to explore the complex mechanism of action of this highly infectious virus. The improved understanding of the structural features of SARS-COV2 (as provided in the supplementary videos) can help to design targeted therapies. Applied mathematical models can help to explore the dynamics of these interactions more accurately Al-Utaibi et al. (2021); Yu et al. (2020); Abdel-Salam et al. (2021); Yu et al. (2021). In this manuscript, we have worked on a model that is linked with SARS-COV2 infection mechanism. The mathematical approach demonstrates how SARS-COV2 is more efficient and adapted to human cells. The model is developed with the aid of the cell-protein interaction studies, available in the literature. The concept of delay has not only proved to be an important, but a deadly weapon of this virus. Our mathematical model features this line of action of SARS-COV2 more accurately.

The rest of the manuscript is organized as follows: In Sect. 2, the mathematical model with delay is presented.

In Sect. 3, the stability analysis, HOPF bifurcation, hybrid genetic algorithm and important cases are presented.

In Sect. 4, important results are discussed and at the end, useful conclusions are drawn.

Model development

To develop the model, we need to first synchronize the biological phenomena with the mathematical modeling procedure. For this the details are provided in videos S1 and S2. The schematic can also be understood with the aid of Fig. 1.

Fig. 1

Schematic depiction of the model

Scheme of the viral adaptation of SARS-CoV2 in the improvement of the infectious process. (A) Acquisition of the adaptive mutation D614G; (B) Structural adaptive modification for the furin cleavage site; (C) Protein S structured from the complex of adaptive modifications that improves the rate of viral infection.

Computational tools have always helped to explore the biological phenomena in a cost effective manner Nutini and Sohail (2020). These methods includes modeling, simulation and forecasting tools. The viral pathology can be interpreted with the aid of mathematical models Belz et al. (2002); Iftikhar et al. (2020). The treatment strategies can be explored with the aid of the mathematical models in an efficient manner, at different scales. Dynamics at cellular, subcelluar and molecular scales can be modeled with the aid of the hybrid modeling approaches Iftikhar et al. (2020).

In this manuscript, we have developed a model, inspired by the work of Pawelek et al. (2012); Sohail and Nutini (2020); Brooke and Prischi (2020); Bhowmik et al. (2020, 2020); Hoffmann et al. (2020); Caufield et al. (2018); Chen et al. (2020); Kleine-Weber et al. (2018) and the references therein. The model is based on the virus concentration, the target cells, the infected cells at levels 1 & 2 (two levels of action are discussed in the introduction (see Fig. 1 as virus infected cells and virus spreading cells) and the IFN signaling proteins.With initial conditionsDescription of variables and parameters is in Tables (1) and (2) and the dynamics can be well understood with the aid of Fig. 2.

Table 1

Description of Compartments

Symbols	Description
V(t)	Virus load
X(t)	Uninfected target cells
Y(t)	Populations of infected cells at first stage
Z(t)	Populations of infected cells at second stage
F(t)	The effect of interferon (IFN)

Table 2

Description of parameter

Symbols	Description
η	Constant infectivity rate of interaction of V(t) with X(t)
a	Transition rate
c_1	Rate of effectiveness in transition
r	Constant rate of F is secreted by Z(t)
b	Virus production rate
γ	Constant degrades rate
κ	Rate of Virus cleared from the cells
ψ	Death rate of infected cells
c_2	Rate of effectiveness in virus production

Fig. 2

Schematic depiction of the model

Positivity of solutions specifies the existence of cells.

Theorem 2.1

Assume that initial solution , , , and , then the solution of model (1) are non-negative .

Proof: From 2nd equation of model all the parameters has positive values, as By integrating we obtainthe above expression shows that X(t) depends on X(0). Therefore, X(t) is positive if X(0) is non negative.

First equation of model (1)Similarly for third, four and fifth equation of the model we have following resultsFrom equation 4, 5, and 6 it is easily seen that if initial solution is non negative then the solution for all values of time t is non negative.

Equilibrium points

The model (1) has infection free equilibrium point, gain by putting right hand side of equation of the model (1) equal to zeroThe linear stability of model is established with method of next-generation operator on model. The reproduction number of model, indicated by , can calculated asThe capability of virus to produce infection or to be uninfected can be analyzed by basic reproductive number . With refers to a decrease in virus production of infected cells where as infection produce due to the increase in virus infected cells production.

Existence of equilibrium points

Theorem 2.2

The model has exclusive endemic equilibrium point if and only if .

Proof: By calculating endemic equilibrium point, we getwhere:

Results

Stability analysis and the Hopf bifurcation

Here we examine qualitative behavior of the model (1) by analyzing local stability of equilibrium points and Hopf bifurcation, which presents the behavior of model (1) by a small change of the solutions as reaction to changes in the particular parameter. As time delays have the significant effect in complexity and dynamics of this model (1), we will assume them as the parameter of bifurcation. Now we examine stability at endemic equilibrium point, the Jacobian matrix at iswhereThe characteristic equation at endemic equilibrium point iswhereThe coefficients areHere, we discuss stability of endemic equilibrium and Hopf bifurcation conditions of the threshold parameters such as and by assuming different cases.

Case 1. When both delay ,and are zero equation (13) becomeEndemic equilibrium is asymptotically stable by Routh-Hurwitz Criteria ifholds, then all the roots are negative. Where and .

Case 2. For and is a real positive number, equation (13) turn out to beWe suppose that there exists real positive number for some value of in such a way that is the root of (18), then we have two equationsAfter simplifying these equation we havewhere the constants areBy rule of signs of Descartes, equation (19) has as a minimum one positive root if and holds.

By eliminating form equation (19) we havewhereDifferentiating equation (18) with respect to delay with the assumption of , then transversality form is obtainwhereThe hopf bifurcation arise for delay if . The above analysis is summarized in following theorem.

Theorem 3.1

Assume that and holds, where delay , in that case, there exist such that is locally asymptotically stable for and unstable for , where in equation (22). Furthermore, at the model (1) undergoes Hopf bifurcation at endemic equilibrium point.

Case 3. When and , in same procedure of case (2), we reach at subsequent theorem.

Theorem 3.2

For model (1) where , in that case, there exist such that is locally asymptotically stable for and unstable for , where in equation (26). Furthermore, at the model (1) undergoes Hopf bifurcation at endemic equilibrium point,where

Case 4. When both and are positive. Then, suppose that as variable and is fixed parameter on stable interval. Assume that there exist a number such that is the root of (13), we obtainAfter simplifying we haveWhere:By applying rule of signs of Descartes equation (29) has minimum one positive root if and holds. we havewithFor Hopf bifurcation will be fixed and differentiate with respect to in equation (28) by putting at ,whereFrom equation (33) if , then Hopf bifurcation occur at .

Theorem 3.3

If and holds with then, there exists such that endemic equilibrium point is asymptotically stable for and , where in (31). Furthermore, the model (1) undergoes Hopf bifurcation at .

Theorem 3.4

If endemic equilibrium point for then, there exists such that endemic equilibrium point is asymptotically stable for and , where in (35). Furthermore, the model (1) undergoes Hopf bifurcation at .

where

Parametric evaluation with hybrid genetic algorithm

A hybrid genetic algorithm combines the power of the genetic algorithm (GA) with the speed of a local optimizer.

The parametric approximation is the most challenging task after designing a mathematical model and after finding the intervals of stability, i.e. the parameters that satisfy the stability criteria. Optimizing parametric values for mathematical models has always remained a great challenge Abdel-Salam et al. (2021).

With the advancement in the field of artificial intelligence and data sciences, the parametric approximation is made easier, keeping in view the stochastic, probabilistic and/or the randomized nature of the real data sets.

In this manuscript, we have used a hybrid optimization tool, partially based on the genetic algorithm, that works for several populations of the parametric mutated genes (sets of values). Matlab platform was utilized for this purpose. Furthermore, the parametric values are selected by keeping in view the intervals imposed by the biological characteristics of the viral process of infection.

A continuous genetic algorithm, that can easily hybridize with the local optimizer, is used during this research. In simple words, the improved values from the genetic algorithm are carried forward by the local optimizer to reduce the computational complexity.

Numerical simulations

We have run some numerical experiments for the understanding of virus control and on the other hand, the bifurcation, linked with the delay.

Figure 3 depicts the role of important parameter b, in understanding the virus spread. For different values of b, we have obtained different dynamics and since the virus replication rate is directly proportional to b, for increased values of b, the virus spread increases and the phase space provides a better understanding of increase in infection, relative to virus load, target cells and the Furin action (see arrow indicating the peak in amplitude). Similarly, Fig. 4 provides information about the change in parameter, linked with the different infection stages (i.e. moving from the compartment of infected cells at first stage to infected cells at second stage). The change in angle of the phase portrait provides useful information about the dynamics.

Fig. 3

Impact of virus reproduction on: a virus load, b phase plot for healthy cells, virus load and Furin

Fig. 4

Impact of infection stages on: a virus load, b phase plot for healthy cells, virus load and Furin

Figure 5 provides useful information about the impact of delay in transmission from one compartment to another, on the virus replication, infected cells and Furin. We can see that for increased delay, as anticipated analytically, there is bifurcation.

Fig. 5

Impact of delay on: a virus load, b phase plot for healthy cells, virus load and Furin

Summary of results

A mathematical model is analyzed with non-negativity of solution, equilibrium points and stability analysis.

Theorem 2.1 shows that the values of compartments is always positive as the parameter is positive.

The Basic reproductive number is obtained. It is calculated by the model of ordinary differential equations, using analytical approach and Matcont numerical approach.

If basic reproductive number , the infection free equilibrium point is stable and infection is completely vanished.

If the basic reproduction number , endemic equilibrium point is stable in feasible interval.

Here we use time delay as parameter of bifurcation to examine Hopf bifurcation.

The non negative endemic equilibrium point is stable when the time delay is very small as time delay increases, the instability occurs that is in accordance with the hopf bifurcation criteria.

Hopf bifurcation is use to find out the instability region in the neighborhood of endemic equilibrium point.

Considering both the D614G mutation and the facilitated action of furin in this process, we assume parameters inclusive of these characteristics.

Discussion

The impact of the SARS-CoV2 virus is devastating mainly due to its speed of infection. The proposed model analyzes: During this research, it is observed that the model is sensitive to the parameters. These parameters were taken from the literature as mentioned in the introduction and the mathematical modeling section. The parameters are responsible for the furin action and SARS-COV2 action. This fact is demonstrated well, with the aid of the numerical simulations, emerging from the Matcont software and genetic algorithm toolbox. The software has the facility for the parametric approximation as well as for the simulations with parametric sweep. The numerical experiments for different values of the parameters and the delay variable are presented in the previous section.

Action of enzyme “Furin” in the speed and spread of the virus.

Presence of D614G mutation (video S1).

Limiting value for , i.e. the interaction rates.

Realistic connection of delay in time with the host and virus interactions.

Importance of delay in the interacting populations of infected cells at second stage, F(t) and the effect of interferon (IFN).

Conclusion

The impact of the CoViD-19 pandemic is devastating and scientific research seek to understand the mechanisms of infection, to create an appropriate vaccine. This paper analyzes the characteristics of the SARS-CoV2 viral infection that shows a fundamental adaptation in the infection process. Arg-Arg-Ala-Arg (RRAR) cleavage site “RRAR” is a cleavage site for the“convertase furin” pro-protein, and is found in the spike protein (S), exclusively in SARS-COV2 virus and is involved in the activation of S protein. In this manuscript, the action of Furin is demonstrated in detail with the aid of the IFN, virus and human cell interaction dynamics. Variable delay helped to link the model with the real dynamics. We conclude that the modeling approach can be further improved by linking it with the forthcoming results from the clinical trials.