The modified Kies-Fréchet distribution: Properties, inference and application

Abbreviations: CAR-T cells: Chimeric antigen receptor T cells; OVs: Oncolytic viruses; TME: Tumor microenvironment; ODE: Ordinary differential equations; eFAST: Extended fourier amplitude sensitivity test

1.   Introduction

An emerging and promising cancer treatment strategy involves the use of chimeric antigen receptor (CAR)-T cells. In this therapy, T cells are extracted from a patient, genetically engineered to express chimeric antigen receptors (CARs), which aid T cells to recognize tumor-derived antigens or neo-antigens and eliminate tumor cells ^[1,2,3], and are infused back to the patient. CAR-T cell therapy has received the Food and Drug Administration (FDA) approval for treating B-cell malignancies ^[4]. Despite the remarkable therapeutic success of CAR-T cell therapy for blood cancers ^[5] and the promising results in treating hematological cancers and CD$19^{+}$ B cell malignancies ^[6,7,8], there are still several challenges that hinder a broader clinical application ^[6]. Some major limitations include reduced identification of tumor-spefic antigens ^[9], immunosuppressive tumor microenvironment (TME) ^[10], and antigen loss in solid tumors ^[11]. Moreover, in some cases, CAR-T cells express certain T cell receptors (TCRs) which are found on normal T cell surfaces ^[12]. This often leads to graft versus host reactions that results in destruction of normal cells. To overcome this challenge, the patient's own T cells are used to generate CAR-T cells. The extraction of T cells from the patient is, however, very costly and sometimes unsafe ^[12,13]. Downregulation or complete loss of target antigens on tumor cell surfaces (also known as antigen escape) is another hurdle, which often leads to tumor escape from T-cell mediated immunity ^[11,13,14], faced by CAR-T cells ^[15,16]. Hence, there is a need for more research to explore the potential of combination therapy that can synergize with CAR-T cell therapy to surmount some of these challenges experienced by CAR-T cells.

Therapy with oncolytic viruses (OVs), which are naturally occurring or genetically engineered viruses that selectively replicate within tumor cells without damage to the normal tissues ^[17,18,19], recently emerged as a promising novel cancer modality that may be further combined with existing therapeutics to enhance their effects. In modern gene-based therapies, OVs have emerged as promising antitumor therapeutic modalities in scientific research and clinical practice ^[20]. In addition to their direct tumor-selective oncolysis, OVs posses many attractive properties that enable them to be used in combination with other treatment modalities ^[12,21,22]. OV can be genetically engineered to specifically lyse tumor cells, while sparing normal healthy cells, by either deleting or inserting new therapeutic transgenes ^[23]. Moreover, OVs can directly or indirectly (or can be genetically modified to) stimulate the immune system ^[24], which can further provide additional antitumor activity. Despite their varying efficacies in treating several cancers, accumulating evidence indicates that OVs been shown to provide a robust antitumoral effect in numerous preclinical studies ^[20,25]. However, these promising results, when used as a standalone treatment, the therapeutic benefits of OVs have so far been limited. In particular, systemic delivery of oncolytic virus has been restricted by antiviral immune response, neutralizing antibodies, and complement proteins ^[26], to mention a few. Even though tumors develop numerous defective or immune escape mechanisms, these mechanisms render them more susceptible to viral infections ^[20]. Furthermore, and most remarkably, OVs are known to heat up "cold" tumors into "hot" tumors, which are easily recognized by adaptive immune cells, especially T cells. Hence, OVs, by their ability to revert the immunosuppressive tumors, to specifically infect and replicate in tumor cells, as well as their capability to activate and recruit immune cells into the TME, are attractive agents for cancer immunotherapy.

Clinical evidence suggests that using OVs offers a great potential in overcoming challenges, such as tumor-mediated immunosuppression, which are faced by T cells or CAR-T cell therapy ^[27]. Additionally, OVs can foster infected tumor cells to express target tumor antigens that can be recognised by CAR-T ^[12]. Furthermore, some OVs have the ability to stimulate CAR-T cells's recruitment into the TME and enhance their ability to kill tumor cells ^[28]. Several other studies established that a combination therapy of OV and CAR-T cells, is a promising treatment modality for treating solid and metastatic tumors ^{[1,2,3,9,11,22,29,30,31]}. However, these studies also highlighted that even though it is quite beneficial to use OVs and CAR-T cells in combination, it is extremely important to find an appropriate time between OVs and CAR-T cells injections. An experimental study indicated that there was a modest tumor reduction when OVs and CAR-T cells were administered at the same time compared to when administered sequentially ^[32]. Most importantly, compelling evidence indicates that timing of OVs and CAR-T cells dosage is an important factor that can impact the treatment outcome ^[33,34,35], even though the strength and timing of OV-induced CAR-T cell recruitment varies depending on virus biology ^[36,37]. In particular, an effective combination therapy, requires that OVs are injected first to lyse tumor cells or modulate immunosuppressive tumor TME for CAR-T cells function ^{[11,33,38,39]}. It is therefore natural to ask whether a single dose or multiple sequential doses of OVs during treatment with CAR-T cells can have any synergistic effect on tumor cell reduction. Of course, this depends on (a) the virus oncolytic potency (viral clearance, infection rate and burst size) and (b) the intensity of the virus-induced synergy (virus-induced increased levels in the migration of CAR-T cells to the TME, and in their killing ability of tumor cells). Thus, the range of OV infection, replication and lytic OV potencies during the combination treatment with CAR-T cells also require further investigation ^[33].

Mathematical models are useful tools to explore these important aspects, they have provided valuable insights in various tumor immunotherapy under preclinical, translational and clinical settings (e.g., see the reviews ^[40,41,42]) and in CAR-T cell therapy (e.g., see the review ^[43]). In recent years, several mathematical models that explore different aspects of tumor growth and progression under CAR-T cell therapy have been developed ^{[44,45,46,47,48,49]}. In particular, a recent work by Le$\acute{o}$n-Triana et al. ^[45] investigated the effects of CAR-T cells targeting two antigens (CD$19$ and tumour-associated antigens) against glioblastoma. Of note, with model simulations, it was shown that injecting large numbers of activated CAR-T cells can overcome the challenge of immunosuppressive TME. The complex dynamic interactions between tumor cells, various immune cells, and OVs have been extensively studied in numerous models ^{[50,51,52,53,54,55]}. To the best of our knowledge, there is only one mathematical model that considered the dynamic effects of a combination therapy of CAR-T cells and OVs. In particular, Walker et al. ^[56] developed two different mathematical models that describe tumor response to cytotoxic activities of CAR-T cells and OVs. With model simulations, the authors illustrated that there are some synergistic effects between CAR-T and OVs with respect to total tumor reduction ^[56]. The model further predicted that the combination therapy of OVs and CAR-T cells resulted in a significant reduction of tumor cells. Moreover, model simulations showed that when co-injected, OVs and CAR-T cells could not completely eliminate tumor cells ^[56]. However, these models ignored the OVs-induced synergistic effect on CAR-T cells' recruitment and their ability to kill that is reported in ^[28]. We note that this synergy can cause CAR-T cells to significantly reduce the infected tumor cell population, which may hinder the replication of OVs. It is therefore of paramount importance to investigate the effect of the reported synergy between OVs and CAR-T cells on the outcome of the combined CAR-T/OV therapy.

For this purpose, we develop a simple mathematical model that describes the local interactions between effector CAR-T cells, tumor cells, and OVs. For the sake of keeping our model simple, our model specifically mimics and considers brief encounter of tumor cells with OVs and CAR-T cells under the treatment scenario of immunodeficient "mouse" model. In particular, we explore whether there is any potential synergism in administrating OVs prior to CAR-T cell therapy. While each treatment modality can still cytoreduce tumor growth by its own, we hypothesize that the efficacy of CAR-T cell therapy can be greatly enhanced by prior administration of OVs. Hence, we shall investigate various dose regimens that may provide robust virus-induced synergistic effects on CAR-T cell cytotoxicity against targeted tumor cells.

2.   Materials and methods

2.1.   Mathematical model

As a first step in understanding the complex dynamic effects of CAR-T cells in oncolytic virotherapy, we devise an ordinary differential equations (ODE) model that focuses mainly on the factors regulating tumor response to a combination therapy of OVs and CAR-T cells. While there exists a wide variety of mathematical models in oncolytic virotherapy, based on ODEs ^{[50,57,58,59]}, partial differential equations (PDEs) ^[60], stochastic differential equations (SDE) ^[61], as well as models based on hybrid computational methods such as multi-scale ^[62] and agent-based modeling approaches ^[63], in the present study, we considered a simple ODE modeling framework with $4$ variables to understand the basic interaction dynamics between tumor cells, CAR-T cells and OVs within the TME. We adopted this simplistic modeling framework because ODE models are, generally speaking, easy to analyze and the resulting solutions are easily tractable ^{[64,65,66,67]}. Most importantly, some significant insights gained from such simple ODE-based models are easily validated with experimental data ^{[18,58,68,69]}, even with clinical data ^[70], and could also provide vital information pertaining to clinical implementation of OV-­based treatment protocols ^[71,72]. The model consists of three cell populations: the uninfected tumor cells ($T_{u}$), tumor cells infected with OVs ($T_{i}$) and effector CAR-T cells ($C_{T}$), and an OV population ($V$). A schematic representation of the local interactions considered in the model is depicted in Figure 1, and the biological description of the model parameters and their values is provided in Table 1.

A system of equations describing the cellular and OV dynamics of the populations considered in this study is given by

In Eq (2.1a), the first term accounts for the logistic growth of tumor cells in the absence of therapy. The second term denotes the infection of tumor cells by OV at the rate $\beta$. It is important to note that in our modeling framework, we have assumed that there are no multiple infections. That is, once a virus has successfully infected a certain tumor cell, then that cell cannot be infected by another virus. Hence, the probability of a virus infecting a tumor solely depends on the number of available uninfected tumor cells. Mathematically, this dependency is modeled through the frequency-dependent kinetics, $\beta T_{u}V/\left(T_{u}+T_{i}\right)$, as done in previous models ^[73,74,75]. Moreover, the choice of this frequency-dependent infection for our model is justified by our assumption that we are modeling interactions between tumor cells, OV and CAR-T cells within the TME. Thus, the OV infection rate is governed mainly by the virus concentration within the TME and the frequency of available uninfected tumor cells. The third term accounts for tumor cell killing by the effector CAR-T cells. The presence of the term $\left(1+\epsilon V\right)$ is to reflect the OV-induced increase in CAR-T's killing of tumor cells. OV infection causes a cascade of subcellular inflammatory events that often result in recruitment of CAR-T cells to the tumor, and enable tumor to be easily infiltrated by effector CAR-T cells ^[33,36].

The increment of infected tumor cell population is represented by the first term in Eq (2.1b). The second term denotes the killing rate of tumor cells by the effector CAR-T cells. The third rate represents the lysis of infected tumor cells with the rate $l_v$.

In Eq (2.1c), the first term, $\lambda_{V}(t)$ (discussed in detail later), denotes the injection of OVs. The second term accounts for an increase in OV concentration in the TME due to the lysis of infected tumor cells with the burst size of $b$. The killing of OVs by effector CAR-T cells is denoted by the third term. The last term represent the viral clearance.

In Eq (2.1d), an administration of adaptive CAR-T cell therapy is represented by the first term, $\lambda_{C_T}(t)$ (discussed in detail later). The second term represents an OV-induced recruitment of CAR-T cells as in ^[28]. The third term denotes that effector CAR-T cells are stimulated by tumor cells during their interactions, which result in proliferation of effector CAR-T cells with constant rate $\rho_{T}$. The fourth term, which includes the natural death, describes inactivation of effector CAR-T cells by tumor cells following their interactions at the rate $\alpha_{T}$.

One disadvantage of many preclinical models is a frequent dosage of OV injections in order to achieve a statistically significant treatment outcome ^[76], as mathematically highlighted in ^[50]. In the present study, we propose a single dosing schedule for each of OV and CAR-T cell based therapies. We propose to use the following hyperbolic tangent type function

It is clear that this function is smooth enough to ensure the well posedness of our model. More importantly, this function has two interesting properties: (a) it tends to $A$ (respectively $B$) for small (respectively large) values of $t$, and (b) it moves from $A$ to $B$ around the time $t = \tau$, with a gradient that is controlled by $k$.

Using this function, we are now able to construct the injection rates $\lambda_V(t)$ (in units: $plaque-forming units (PFU)\cdot cell^{-1} \cdot day^{-1}$), in Eq (2.1c), and $\lambda_{C_T}(t)$ (in units: $cells\cdot day^{-1}$), in Eq (2.1d), as follows:

The functions $\lambda_V(t), \ \lambda_{C_{T}}(t)$, start and stay at the value $\lambda_{0_X}$ until time $t$ is close to $\tau_{_1X}$ where they move to $\lambda_{0_X}$ and remains at this value until around $t = \tau_{2_X}$ where they move to $0$.

To have a profile that a is typical of an OV or CAR-T single dose injections, we chose $\lambda_{0_V} = 0\ PFU\cdot cell^{-1} \cdot day^{-1}, \ \lambda_{0_{C_{T}}} = 0\ cells \cdot day^{-1}, \ \lambda_{m_V} = 10^{6}\ PFU\cdot cell^{-1} \cdot day^{-1}, \ k_{1_X} = k_{2_X} = 1.5\ unitless, \ \lambda_{m_{C_{T}}} = 2\times 10^{6}\ cells \cdot day^{-1}$ and $\tau_{1_V} = 40\ day, \ \tau_{2_V} = 42\ day$ and $\tau_{2_{C_{T}}} = \left(\tau_{1_{C_{T}}}+2\right)\ day$, where $\tau_{1_{C_{T}}}$ takes various values that accommodate our simulation scenarios. The resulting profile of a single dose injection of OV is given in Figure 2.

For the sake of our simulations, we assume that at time $t = 0$ the tumor is at its maximum carrying capacity, with no infected cells, viruses or CAR-T cells. This ensures that the tumor size is always the same regardless of the times of the initial OV and CAR-T doses are. Note also that the maximum tumor carrying capacity, $K = 1.404\times 10^{8}$ cells, considered in this study, is below the carrying capacity of $K = 3 \times 10^9$, which is equivalent to a tumor size of $3000\ mm^3$ that is lethal to mice ^[73]. Accordingly, the initial conditions of the system (Eq (2.1a)–(2.1d)) are

We note that when $\epsilon = \xi_V = 0$, model (2.1) describes the dynamics of tumor cells under the combined CAR-T/OV therapy without synergy. Compared to the ODE model in ^[56], our model without virus-induced synergy, includes the OV terms and excludes CAR-T cells circulating in the blood, and the overall inflammatory response. We note that there are other plausible possible formulations to mathematically characterize virus-induced synergy between OV and CAR-T cells, but we choose this simple formulation as there is no appropriate experimental data to validate our model about the kinetics of the parameter $\epsilon$.

2.2.   Parameter values

Since our modeling framework is based on immunodeficient mouse models, for CAR-T cell therapy, the baseline parameter values in Table 1 were taken directly from CAR-T cell immunotherapy models in ^[44,45,77]. For oncolytic virotherapy, the baseline parameter values were taken from various mathematical models ^{[51,66,71,78,79,80]} because not all parameter values were available in one study. Indeed, various oncolytic viruses result in different treatment outcomes, partly because of their intrinsic or genetically engineered characteristics. For example, considering a vesicular stomatitis virus (VSV), infected tumor lysis rate $(l_{v})$ is $11$ per day in ^[71], $1$ per day in ^[79], while $0.475$ per day is reported in ^[80].

Thus, there is a diverse, yet very scarce information regarding lysis rate of infected tumor cells by OV. Similarly, it is very difficult to find an appropriate value for viral particles than can trigger infection within the tumor. Thus, our model has a few "free" parameters which we vary within a plausible range to assess their impart on model predictions. These free parameters are infection rate $\beta$, lysis rate of infected tumor cells $l_{v}$, rate of CAR-T cell recruitment by OVs $\xi_{V}$ and viral burst size $b$.

3.   Results

The aim of this work is to investigate how tumor cells respond to a combination therapy of CAR-T cells and OVs, and provide a quantitative description of local effects of CAR-T cells during oncolytic virotherapy. Especially, we explore the effect of the virus oncolytic potency and the intensity of the virus-induced synergy on the output of the combined therapy. In particular, we investigate the effects of (ⅰ) infection rate, (ⅱ) burst size, (ⅲ) viral clearance, (ⅳ) virus-induced increased levels in the migration of CAR-T cells to the TME, and (ⅳ) virus-induced increased levels in CART-T cells' killing ability of tumor cells.

We first hypothesize that timing of OV injections is a major determinant of a combination therapy with CAR-T cells. A major reason for timely injections of OVs and CAR-T cells is to prevent premature OV clearance by effector CAR-T cells. For effective combination therapy, it is also equally crucial to minimize tumor escape by injecting CAR-T cells at an appropriate time to provide complementary anti-tumor effect. Thus, the results in this study are categorized into two major components: Model simulations for combination therapy without virus-induced synergistic effects and Model simulations for combination therapy with virus-induced synergistic effects. The numerical simulations were performed in Matlab (R2017b) using ode23s solver. We now discuss these results in detail.

3.1.   Global sensitivity analysis

We performed a global sensitivity analysis using the extended Fourier Amplitude Sensitivity Test (eFAST) ^[84], a variance-based sensitivity method, to understand how the uncertainty in parameter values impacts the model predictions with respect to the total tumor cell population. The parameter values were varied from half to twice the baseline values in Table 1 to account for a possible range of uncertainty in the model parameter values. eFAST provides two sensitivity indices, namely the first-order indices, $S_{i}$, and the total-order indices, $S_{ti}$. The first-order indices ($S_{i}$) measure the influence of a particular parameter on the variance of the total tumor cell population (the model output) without the local interactions with other model parameters, while the total-order indices ($S_{i}$) measure the influence of a particular parameter on the total tumor cell population with possible interactions with other model parameters. Given these two indices, it is easy to quantify whether the total tumor reduction occurs due to a direct influence of the given parameter or through the influence of other model parameters. To compute the sensitivity indices with the eFAST method, we used a sample size of $4000$ from the normal distribution of each parameter to examine the uncertainty and sensitivity of our model predictions to variations in model parameter values. We implemented the eFAST method using the generic MATLAB code developed by Kirschner et al. ^[110], and the results are given in Figure 3. We find that the reductions in the total tumor cell population are also heavily affected by the virus infection rate ($\beta$) and the rate at which infected tumor cells are lyzed by the OVs ($l_{v}$) at day $40$. The tumor growth rate ($r$) also significantly affects the model predictions (Figure 3(A)), at day $10$ and $40$. Intuitively, this means that tumor aggressiveness (governed by the proliferation rate $r$) plays an important role in determining the likelihood of the therapy success or failure ^{[50,68,85,86]}. On the other hand, the virus-induced synergistic rate ($\epsilon$) and the virus-induced CAR-T cell recruitment rate ($\xi_{V}$), parameters which we have assumed their values in our model, have very little influence on the model predictions, since they have lower $S_{i}$ and $S_{ti}$ values than the dummy variable (last column bar) in Figure 3(A). Therefore, the uncertainty of these two parameters does not significantly affect the robustness of the predictions. Since the carrying capacity, $K$, was the most influential parameter in the model, as evidenced by both the first-order indices and total-order indices, in the numerical simulations, we set the tumor size to its maximum carrying capacity at time $t = 0$, with no infected cells, viruses or CAR-T cells. This assumption ensures that the total tumor size is the same regardless of the initial CAR-T and OV doses, and hence would allow us to quantify the effects of most sensitive parameters and the order of dosages of the two therapies. Collectively, the sensitivity analysis illustrates that, in addition to the tumor carrying capacity, the model is mostly affected by few virus related parameters. Thus, we explicitly study the influence of such parameters on total tumor cell numbers by individually varying these parameters to predict their effects on the combination therapy.

3.2.   Combination therapy without virus-induced synergistic effects

It is argued that CAR T-cells are less likely to synergistically benefit from OV therapy if OV infection induces an apoptotic tumor cell death, rather than immunogenic cell death ^[31]. We investigated this proposition by investigating the behavior of model, Eq (2.1a)–(2.1d), without virus-induced synergy, i.e., $\epsilon = \xi_V = 0$. This therapeutic scenario resembles the case where OVs induce apoptotic cell death (i.e., tumor cell death does not induce any T cell responses). Here, in order to mimic a more plausible dosing schedule, we propose a dosing schedule that allows for sufficient virus propagation in infected cells, depicted in Figure 4.

Note that in Figure 4, the dosing day $42$ for CAR-T cells is the same time used in the experiments in ^[77], while day $17$ is an adhoc time for OV injection. Indeed, for the combination therapy with adoptive T cells, OVs can be injected as early as day $7$ post tumor implantation in mice ^[87].

3.2.1.   OVs with low infectivity can be effectively complemented with adjuvant CAR-T cell therapy

As suggested in the literature that OVs should be administered first to change the TME to a state that favors CAR-T cell function and permissiveness, and directly induce oncolytic effect on tumor cells ^[33], we explored the effects of varying OV infection strength (described by parameter $\beta$ in Eq (2.1a)) on tumor growth. We believe that the strength of OV infection can determine whether adding an adjuvant therapy of CAR-T cells can lead to improved or deleterious effects on tumor growth. Indeed, OV infection rate plays a vital role in oncolytic virotherapy ^[71,88]. The therapeutic effects of varying OV infection on tumor cell killing are shown in Figure 5.

As can be seen in Figure 5(A) when infection rate is low (i.e., $\beta = 9.0\times 10^{-4}$), OVs alone can ultimately lead to tumor elimination as expected, but it takes unrealistically long time to reach a tumor-free state for given duration of in silico experiment, indicating the need for an adjuvant therapy. It should, however, be noted that if OVs infect tumor cells at a little higher rate (i.e., $\beta = 6.6\times 10^{-4}$), then there is a notable and observable difference in the number of infected tumor cells (Figure 5(D)) compared to those receiving slower infecting OVs (in Figure 5(C)), indicating that change in the infection rate seems to bring significant effect on tumor growth. In all cases, we note that CAR-T cells leads to rapid elimination of tumor cells.

Similarly, it is equally important to examine the effects of varying infection rates on tumor cell growth. If OVs have a lower infection rate, then many free OVs would remain vulnerable to CAR-T cell clearance, as depicted in Figure 6. Of note, in Figure 6(A), (B), there is a relatively slow increasing virus accumulation in the TME, which is swiftly interrupted by injection of effector CAR-T cells. While injection of CAR-T cells on day $42$ (see Figure 6(C), (D)) seems not to provide a notable clearance of free OVs within the TME, we note that their cytotoxic effect is more prudent on tumor cells (see Figure 5(A), (B)).

More importantly, consistent with experimental studies ^[77], we note that the effector CAR-T cells rapidly decay to zero in the absence of OVs and tumor cells. We should emphasize that decaying to zero means that either the CAR-T cells are eliminated from the TME due to natural death or they are becoming memory CAR-T cells, and would be readily available upon tumor re-challenge, as shown in ^[44,89].

3.2.2.   A single CAR-T cell dosage is inadequate to eliminate tumor

As a first step, we investigated where a single dose of $2\times 10^{6}$ cells is sufficient to drive the total tumor cell population to extinction. Figure 7 indicates that a single dose of $2\times 10^{6}$ CAR-T cells sufficed to reduce the total tumor load to a smaller size in a few days; however, the tumor rapidly regrew to its maximum carrying capacity. This implies that a single dose of CAR-T cells can lead only to a small transient tumor reduction, and thus the therapy cannot control tumor growth in this in silico experiment (see Figure 7). We should emphasize that this observation could be related to the particular choice of parameter values and initial conditions used in our model, and would significantly change under different parameter values and initial conditions.

Thus, an initial dose of CAR-T cells results in small reduction in the total tumor cell population (Figure 7), necessitating the need for additional CAR-T cell doses to control tumor growth in silico. It imporatnt to note that this observation is related to the particular choice of parameter values used in this study, and the outcomes are likely change under different conditions.

Adding a single dose of OV, we observe in Figure 8 that using a single dose of CAR-T cells given at any of the times given in the legend, the total tumor size converges to zero.

Interestingly, in Figure 8, we observe that using a single dose of OV and single dose of CAR-T cells given at any of the times given in the legend, the total tumor size reaches to zero.

In summary, using our proposed model, without virus-induced synergism, we explored diverse scenarios which can lead to eradication of tumor cell populations. In particular, it was shown that variations in OV infection rates may determine if CAR-T cells can augment or inhibit virus oncolysis. Depending on the dosage and timing of OV administration, if the virus has lower infection rate, then its propagation within the tumor cell or TME is likely to be halted by effector CAR-T cells. On the other hand, if the virus has higher infection rates, then it is likely to eradicate all tumor cells prior to injection of CAR-T cells. Collectively, these results illustrate that combating tumor burden with a combination therapy of OVs and CAR-T cells does not always provide superior therapeutic outcome as CAR-T cells can lead to either augmentative or deleterious effects. Furthermore, we showed that a single dose of CAR-T cells is not enough to eradicate tumor cells.

3.3.   Combination therapy with virus-induced synergistic effect

We now use our model, Eq (2.1a)–(2.1d), to investigate where OVs can co-operate with CAR-T cells by inducing immunogenic tumor cell death, rather than apoptotic tumor cell death considered previously. In particular, we investigate whether there is any synergism between OV and CAR-T cells in killing tumor cells and how the order of OV and/ or CAR-T cell administration impact the overall combination therapeutic outcome.

3.3.1.   CAR-T cells energize with OVs when virus concentration is high

To determine whether OVs can synergistically enhance CAR-T cell lytic activity against tumor cells, we simulated the treatment scenario where tumor is treated with a fixed dose of each therapy. In particular, we propose a treatment scenario where tumor is treated with one dose of $10^{6})$ PFU on day $40$ and same fixed doses of $2\times 10^{6}$ CAR-T cells at various time points (i.e., before, concurrent with, and after OV administration). This treatment setting is motivated by, but sightly different from, the in vivo experiment done by Krabbe et al. ^[87] (where the antigen-specific OTI T cells are injected once and oncolytic vesicular stomatitis virus (VSV), VSV-NDV, is injected prior and posterior to OTI T cell injection). Figure 9 illustrates examples of the model's (including the synergy parameter, (2.1a)–(2.1d)) simulations where a potential synergistic effect between OVs and CAR-T cell therapies is considered with respect to dosage and timing of each therapy.

As can be seen in Figure 9(A), CAR-T cell monotherapy is unable to control tumor growth, applied in day $10$ and $20$, but when applied in combination with OVs, on day $40$, there is a significant increase in cytotoxicity of tumor cells, consistent with recent experimental results in ^[87]. It is also important to note that because OVs are also susceptible to CAR-T cell clearance, there would be a rapid decline in free virus concentration within the TME (Figure 9(B)). Since the tumor is not completely eliminated, we observe that the same fixed subsequent dosage of CAR-T cells is still unable to reduce tumor progression, and tumor escapes to some lower stable steady state (Figure 9(A)), though below a lethal threshold of $1 \times 10^{10}$ cells considered in ^[44] or $3 \times 10^{9}$ cells for mice ^[73]. In agreement with the findings in ^[90], depending on the initial tumor size, it is plausible that the total tumor size can be reduced to a smaller size than the carrying. This result further suggests that complete tumor elimination can be observed via increase in CAR-T cell dosages, in agreement with numerical simulations in ^[44]. Of particular note, in Figure 9(A), the maximum tumor reduction is observed when virus concentration is high, indicating synergistic therapeutic effects as a result of combination therapy (i.e., co-administration of OVs and CAR-T cells). While this result seems counter-intuitive because one may expect that OV infection would rapidly be forestalled by CAR-T cells ^[31], the observation in Figure 9 is consistent with many in vitro and in vivo studies that consider a combination treatment of OVs and CAR-T cells ^[2,3,11,87]. One possible explanation for this increased antitumor activity is that OV infection results in immunogenic cell death of tumor cells, which often results in recruitment of T cells to the TME ^{[29,31,91,92,93]}. In other words, this means that since not all OVs will be cleared by CAR-T cells, during simultaneous administration, then the resulting death of infected tumor cells, may further recruit CAR-T cells towards to tumor cells. It is also important to note that there exists a non-monotinic regrowth of the total tumor (Figure 9(A)) with CAR-T doses at day $100$ and $130$ (purple, light blue) due to the small increases in the virus populations depicted in Figure 9(B). This observation could be attributed to low viral replication rate (which is below a certain threshold ^[94,95]), as a result of limited virus infection propagation by CAR-T cells.

Notably, we see in Figure 9 that the combined therapy did not succeed in eliminating all tumor cell populations. This could be explained by the fact that the synergistic effect of Ovs on CAR-T cells, enhances CAR-T cell proliferation and their antitumor activity, which may result in killing of infected tumor cells and hence hinder OV oncolysis. This gives a clear guidance why there is always a need for sequencial administration of multiple doses of CAR-T cells and OVs. The dosage and timing of the OVs and CAR-T cells governs the potential synergism that can result from the combination therapy. In particular, the dosing regime where simultaneous administration of OV and CAR-T cells is applied has the greatest overall tumor reduction compared to the dosing regimen where CAR-T cells are administered first or after simultaneous dosage of OVs and CAR-T cells.

Taken together, these results demonstrate not only the potential of simultaneously administrating CAR-T cells and OVs to obtain an improved therapeutic effect, but also highlight the importance of designing new dosing strategies to simultaneously administer CAR-T cells and OVs at different times. Moreover, it seems plausible that during the co-administration of CAR-T cells and OVs, there is some potential virus-induced synergism that enhances CAR-T cell antitumor activity. While this might be a favorable treatment outcome, it should be noted that if the CAR-T cell induction is too strong, then CAR-T cells can potentially inhibit OV infection propagation by eliminating both infected tumor cells and free OVs within the TME, as was shown with the natural killer (NK) cells by Senekal et al.^[66]. Furthermore, these results highlight the importance of mathematical modeling as a tool that can be used to explore plausible therapeutic outcomes prior to in vitro or in vivo studies.

3.3.2.   Investigating the effects of the virus parameters

In this section, we explore the effect of the virus oncolytic potency and the intensity of the virus-induced synergy on the output of the combined therapy. In particular, we investigate the effects of (ⅰ) infection rate, (ⅱ) viral clearance, (ⅲ) burst size, (ⅳ) virus-induced recruitment of CAR-T cells to the TME, and (ⅳ) virus-induced increased levels in CAR-T cells' killing ability of tumor cells.

As suggested in the literature that OVs should be administered first to change the TME to a state that favors CAR-T cell function and permissiveness, and directly induce oncolytic effect on tumor cells ^[33], we first studied the effects of OV infection strength (i.e., as described by the parameter $\beta$ in Eq (2.1a) on tumor growth. We believe that the strength of OV infection can determine whether adding an adjuvant therapy of CAR-T cells can lead to improved or deleterious effects on tumor growth. Indeed, OV infection rate plays a vital role in oncolytic virotherapy ^[71,88]. We examine the effects of strength virus infection rate by varying the infectivity rate $\beta$ and examined the therapeutic impact on tumor cell killing in Figure 10.

In agreement with other recent theoretical models ^[66,96], if OVs infect tumor cells at a little higher rate (i.e., $\beta = 3.7\times 10^{-4}$), then there is a more observable difference in tumor cell reduction (Figure 10(A)) compared to those receiving slower infecting OVs (in Figure 10(B)), indicating that small changes in infection rates of OVs can have significant effects on tumor growth dynamics.

Another significant parameter that often defines the success of the therapy with OVs is the burst size, which is the number of infectious viruses produced by lyzed infected tumor cell. Figure 11 illustrates the effects of varying the burst sizes on total tumor cell population.

As demonstrated in Figure 11(A), and consistant with other biological ^[97] and mathematical findings ^{[50,51,78,96]}, OVs with higher burst sizes have a greater antitumor impact. The rationale for this is that when OV burst size is high, then there will be more free OVs within the TME, and which, consequently, are likely to infect more tumor cells. It is worth mentioning that one vital attribute of OVs is that they can easily be genetically modified to produce higher viral titers ^[98], and our simulation results further highlight the essence of manipulating the virus burst size within the biologically feasible ranges.

Decay rate of OV, $\omega$, with respect to therapy outcomes, is another important parameter worth considering. As demonstrate in Figure 12(A), if the OV has a very low virus death rate $\omega = 0.9$, there is a greater antitumor impact on tumor, compare to higher virus death rate $\omega = 1.2$, which results in minimal reduction of tumor load from equilibrium (Figure 12(B)). The longer the lifespan of an OV, the higher the lytic effect on tumor cell population.

As can be seen in Figure 12(A), even though there is an immediate tumor reduction, upon sequential injections of effector CAR-T cells at day $10$ and $20$, effective tumor control by CAR-T cells is observed after a long time period (around days $150$ (blue curve) and $245$ (yellow curve) in Figure 12(A)). One possible explanation for this occurrence is that upon the joint administration of CAR-T cells and OVs at day $40$, antiviral activity of CAR-T cells was more robust than antitumor activity; thus, taking a longer time to effectively impact the total tumor load at equilibrium. In a nutshell, our study suggests that a significant reduction of the total tumor load (though it occurred on a much longer time scale) is also achievable with early time points administration of CAR-T cell monotherapy.

Oncolytic virus-induced potency is one of the key determinants of successful oncolytic virotherapy ^[50,99,100], thus, in Figure 13, we investigated whether an increase in the rate OV kills tumor cells, $l_{V}$, and the rate of virus-induced CAR-T cell antitumor activity, $\epsilon$, can improve or inhibit the combined antitumor effect of OV and CAR-T cells.

Of particular note, less lytic OVs do result in more tumor reduction (Figure 13(A)), as compared to OVs with high lytic ability that are able to destroy a large portion of tumor cells (Figure 13(B)). On the same hand, and most interestingly, we also observe that when there is high virus-induced synergistic killing of tumor cells by CAR-T cells, there is a greater total tumor cell reduction (Figure 13(C)), compared to low virus-induced synergism (Figure 13(D)). It is important to note that while this effect (of virus-induced synergism) is therapeutically exciting, it is by far less significant compared to the therapy without virus-induced synergistic effect (Figure 8). Taken together, these results suggest that high virus-induced tumor killing by CAR-T cells may not lead to successful combination therapy immunogenic OVs and effector CAR-T cells. We should emphasize that higher virus-induced synergy can lead to premature and rapid antiviral CAR-T cell responses, if the OV immunogenic strength is too strong as shown in ^[66], which in turn would inhibit OV oncolysis. This explains why the combination therapy without virus-induced synergy led to tumor elimination (Figure 8), while the one with virus-induced synergism did not Figure 13(C), (D).

The effect of recruiting effector CAR-R cells to the TME, is investigated in Figure 14.

Consistent with experimental studies ^[98,101,102] and mathematical models ^[66], in Figure 14(B), we note that OVs that immunogenetically recruit immune cells at higher rates do not offer optimal therapeutic benefits, Nonetheless, the effective recruitment of CAR T-cells into the TME dependents mainly on the cascade of events leading to diverse chemokine signalling pathways ^[32].

We end this section by stating the inherent limitation of the present study: Current limitations of the present study include the exclusion of tumor antigen heterogeneity when OVs and CAR-T cells are administered and fixed dosage regimes. Expression of diverse tumor antigens is one of the major hurdles of CAR-T cell based therapy against solid tumors ^[103,104]. While the mathematical model proposed in this study could be utilized to design and generate in silico trails to explore optimal therapy combinations, it did not consider the possibility of fractionated dosing regimens and timings of the CAR-T cells, which might reduce the possibility of CAR-T cell toxicities. Indeed, fractionated dosing regime for CAR-T cells was shown to reduce observable toxicities while still resulting in more effective CAR-T cell treatment ^[105].

4.   Conclusions

While CAR-T cell therapy has recently demonstrated promising results against many hematological malignancies, its clinical application against solid tumor still faces many challenges, including antigen loss or heterogeneity, immunosuppressive TME and suboptimal CAR-T cell trafficking. As an attempt to improve the efficacy of CAR-T therapy, combination therapy of OVs and CAR-T cells has recently proven to be one of the promising solutions for tackling some barriers of CAR-T cell therapy. Although it has recently been shown that OVs can aid CAR-T cells in various ways, many studies are still needed to demonstrate whether an apoptotic or immunogenic cell death induced by OVs can synergize with CAR-T cells to potently kill tumor cells. In this article, we developed a simple compartmental model, with virus-induced synergism, to describe the local interactions between OVs, tumor and effector CAR-T cells with the TME.

First, using the model without synergy parameter, we explored the potential effects of administering OVs first, then followed by CAR-T cells, we noted that the level at which CAR-T cells can augment virus oncolysis depends on the strength of virus virulence (i.e., whether OVs are strongly or weakly lytic) and the dose order (i.e., whether OVs are administered before CAR-T cells administration). Second, by including the synergy parameter in our model, we further investigated whether the is any potential virus-induced synergism by considering a dosing regime where CAR-T cells are administered prior to, concurrent with OVs and posterior to simultaneous OV and CAR-T cells administration. Our sensitivity analysis (based on eFAST) results revealed that the tumor carrying capacity is a most vital parameter (indicated by high sensitivities in both first-order and total-order indices) that can significantly impact the therapy outcomes. From model simulations, it was shown that simultaneous administration of OVs and CAR-T cells considerably decreased tumor burden compared to CAR-T cell monotherapy, applied prior to or posterior to simultaneous injections of OVs and CAR-T cells. Considering how the ordering of the dosages and timings of OVs and CAR-T cells resulted in different tumor growth dynamics, it is advisable to administer OVs first and/or concurrent with CAR-T cells, to prime or sensitize CAR-T cells towards tumor. This sensitization can result in greater antitumoral effect—due to a cascade of intra-cellular signaling events mediated by OV oncolysis which would ultimately recruit CAR-T cells to TME. Moreover, and notably, our simulations illustrate that tumor response to the combination therapy of OVs and CAR-T cells is directly dependent on the intensity of virus-induced synergistic tumor killing by effector CAR-T cells, which, in turn, would decrease the efficacy of the combination treatment. Taken together, our results suggest that virus-induced synergistic effects of CAR-T cell on tumor cell killing during the treatment with OVs can potentially impede treatment success.

Although the results attained in this study are informative on how synergistically OVs can induce the direct effect of CAR T-cell function on tumor growth, there are several limitations and simplifications that suggest follow up studies. The fact that the total tumor cell killing by CAR-T cells is minimal in the present model, while CAR-T cells are known to have much antitumor effect in the literature ^[5,15,16,30], naturally limits the model applicability and suggests further refinement of the model's assumptions and complexity. Limited CAR-T cell effect is related to the killing rate of tumor cells by CAR-T cells, assumed here to have a fixed value of $10^{-7}$ cells per day across all tumor cells, which is in the range obtained experimentally in ^[44]. While this simplification was made in the spirit of keeping our model simple, in patient or mice populations, there could be some level of variation in the killing rates of CAR-T cell due to diverse antigen expression by both uninfected and infected tumor cells ^[106]. Implicitly, in the current model, we assumed that the CAR-T cells are able to target one antigen which is present to both uninfected and infected tumor cells. Recent advances in CAR-T cell designs focus on making CAR-T cells to target multiple antigens ^[33,45,107]. In future work, it would be interesting to examine how variations in the CAR-T cell killing rates affect treatment success. Another important simplification that limits the present study is that our modeling framework does not account for tumor cell heterogeneity and spatial variations in CAR T-cells, or OVs concentration. Due to the diverse nature of heterogeneity of tumor cells, intrinsic or induced by OV infections, spatial distributions and heterogeneity of tumor cells may be important considerations for future work. Furthermore, it would be interesting to explicitly investigate not only the non-localized interactions between tumor cells, CAR-T cells, and OVs, which are regulated by specific chemokine signaling pathways, but also the effects of temporal dependencies of CAR-T cell antitumor activity on tumor growth. Finally, a future direction of the present work would examine a wide range of CAR-T/OVs mechanisms (e.g., the optimal dosage that leads to maximum tumor reduction as done in ^[108,109]) that can lead to optimal virus-induced synergy necessary to mount potent CAR-T cell antitumor response against solid tumors.

Acknowledgments

We greatly thank the anonymous reviewers whose insights, valuable comments and suggestions have improved our manuscript. R. Ouifki acknowledges the support from the DST/NRF SARChI Chair in Mathematical Models and Methods in Biosciences and Bioengineering at the University of Pretoria. L. de Pillis acknowledges the support of the Norman F. Sprague, Jr., Chair of Life Sciences and the Department of Mathematics of Harvey Mudd College.

Conflict of interest

The authors declare there is no conflict of interest.