TREATMENT STRATEGIES FOR COMBINING IMMUNOSTIMULATORY ONCOLYTIC VIRUS THERAPEUTICS WITH DENDRITIC CELL INJECTIONS

. Oncolytic viruses (OVs) are used to treat cancer, as they selectively replicate inside of and lyse tumor cells. The eﬃcacy of this process is limited and new OVs are being designed to mediate tumor cell release of cytokines and co-stimulatory molecules, which attract cytotoxic T cells to target tumor cells, thus increasing the tumor-killing eﬀects of OVs. To further pro-mote treatment eﬃcacy, OVs can be combined with other treatments, such as was done by Huang et al. , who showed that combining OV injections with dendritic cell (DC) injections was a more eﬀective treatment than either treatment alone. To further investigate this combination, we built a mathematical model consisting of a system of ordinary diﬀerential equations and ﬁt the model to the hierarchical data provided from Huang et al. We used the model to determine the eﬀect of varying doses of OV and DC injections and to test alternative treatment strategies. We found that the DC dose given in Huang et al. was near a bifurcation point and that a slightly larger dose could cause complete eradication of the tumor. Further, the model results suggest that it is more eﬀective to treat a tumor with immunostimulatory oncolytic viruses ﬁrst and then follow-up with a sequence of DCs than to alternate OV and DC injections. This protocol, which was et al. , allows the infection to initially thrive before the immune response is enhanced. Taken together, our work shows how the ordering, temporal spacing, and dosage of OV and DC can be chosen to maximize eﬃcacy and to potentially eliminate tumors altogether.

Abstract. Oncolytic viruses (OVs) are used to treat cancer, as they selectively replicate inside of and lyse tumor cells. The efficacy of this process is limited and new OVs are being designed to mediate tumor cell release of cytokines and co-stimulatory molecules, which attract cytotoxic T cells to target tumor cells, thus increasing the tumor-killing effects of OVs. To further promote treatment efficacy, OVs can be combined with other treatments, such as was done by Huang et al., who showed that combining OV injections with dendritic cell (DC) injections was a more effective treatment than either treatment alone. To further investigate this combination, we built a mathematical model consisting of a system of ordinary differential equations and fit the model to the hierarchical data provided from Huang et al. We used the model to determine the effect of varying doses of OV and DC injections and to test alternative treatment strategies. We found that the DC dose given in Huang et al. was near a bifurcation point and that a slightly larger dose could cause complete eradication of the tumor. Further, the model results suggest that it is more effective to treat a tumor with immunostimulatory oncolytic viruses first and then follow-up with a sequence of DCs than to alternate OV and DC injections. This protocol, which was not considered in the experiments of Huang et al., allows the infection to initially thrive before the immune response is enhanced. Taken together, our work shows how the ordering, temporal spacing, and dosage of OV and DC can be chosen to maximize efficacy and to potentially eliminate tumors altogether.
1. Introduction. Certain cancers do not respond or become refractory to traditional systemic therapies such as chemotherapy, requiring the development of novel treatments [1,13,31]. One such treatment is oncolytic virotherapy [1,31,36,46]. Oncolytic viruses (OV) selectively replicate inside of and kill tumor cells while sparing normal tissues. However, viral replication and lysis are necessary but often not sufficient to completely and durably eliminate a tumor. Novel engineering approaches that combine viral oncolysis with production of immunostimulatory molecules, as well as dendritic cell (DC) injections, are under development with the goal of improving the efficacy of treatment [3,18,23,36].
To this end, certain viruses are being genetically engineered to foster expression of pro-inflammatory cytokines and co-stimulatory molecules which, when released, induce a targeted immune attack against the tumor [15,30]. When coupled with viral oncolysis, this approach produces a robust anti-tumor immune effect.
Additionally, the immune system can be bolstered by the injection of DCs that have been found to act synergistically with cytokine-expressing OVs to produce potent antitumor and immune-system boosting effects [18,36]. Pre-clinical experiments by Huang et al. show that combining cytokine-expressing OVs with DCs is a promising new therapeutic approach that deserves more attention [18].
Currently, there is no consensus about optimal choices for cytokine expression or dosing strategies once a particular viral platform is chosen [27,37,47]. Moreover, including DC injections into the therapy protocol makes the dosing problem exponentially more complicated.
Mathematical models can expedite pre-clinical and clinical research by providing insight into the dynamics of these complex biological systems. Differential equations models paired with experimental data can improve the design and administration of experiments and treatments [14,16,24,26,42,43]. To date, model-based development has been successful in many cases including optimized treatment schedules for measles virus [5]; interferon-evasion mechanisms for vesicular stomatitis virus (VSV) [26]; ideal combinations of adenovirus infection and MEK inhibition [3]; viral burst size needed when combined with cyclophosphamide to reduce glioma tumor burden [16]; and personalized prostate cancer vaccine regimens [25].
Besides furthering clinical experiments by providing an inexpensive method for optimizing virotherapy, mathematical models can also reveal the inner workings of the rich dynamics of a complicated biological system. Mathematical models have been used to reveal oscillations in tumor growth [12], to determine how initial conditions affect the efficacy of OV therapy on melanoma [35], and to determine factors that inhibit and enhance OV transmission through a tumor site [32]. Several mathematical models have also been used to investigate the complex dynamics that arise between tumor cells, oncolytic viruses, and the immune system [41,45]. Other recent models have focused on tumor-immune interactions with DC vaccination [9,33,34].
This article explains our formulation of a mathematical model that includes cytokine and co-stimulatory molecule expressing OVs in combination with DC injections. In particular, we model an OV approach developed experimentally by Yun and colleagues, that combines an adenovirus, which mediates the production and release of intereukin-12 (IL-12) and 4-1BB ligand (4-1BBL), with DC injections.
After model development, we explicate the methods used to fit the model to murine tumor volume data from Yun's lab. The model fitting is done hierarchically, starting with basic tumor growth, advancing through addition of the adenovirus, then adding immunostimulatory cytokines and molecules, and ending with the addition of DC injections. Next, we discuss results from simulations including goodness of model fits, variations in dosing strategies, including different strategies for combining injections of OV and DC, and varying dosing periods. Finally, we discuss results from simulations describing how alterations in the virus affects final tumor size. Overall, we show that OV therapy in combination with DC injections can be an effective treatment for cancer.

Model.
2.1. Model equations and assumptions. In this work, we developed a mathematical model of tumor dynamics in response to treatment by oncolytic adenoviruses and DC injections with the goal of investigating different treatment strategies. For this study, we extended a model that we developed in a previous work [22], that we used to study the dynamics of immunostimulatory oncolytic adenoviruses. In addition to immunostimulatory OV therapy, our new model can accommodate DC vaccine therapy. Two of the oncolytic adenoviruses modeled here mediate expression of either 4-1BBL or IL-12 singularly, and the third virus mediates the release of both by infected tumor cells. This molecule, 4-1BBL, and cytokine, IL-12, recruit and stimulate T cells and antigen-presenting cells (APCs), that are part of the immune attack to eradicate tumor cells. Our model is formulated as the following system of ODEs: where U is the number of uninfected tumor cells, I is the number of infected tumor cells, V is the number of virions, T is the number of T cells at the tumor site, A is the population of endogenous APCs in the tumour-draining lymph nodes, and D is the population of autologous injected DCs. In addition, N = U + I + T is the total population of tumor and T cells at the tumor site. A diagram of the model is given in Figure 1. The mechanisms modeled by (1)-(6) are described as follows: • In (1), uninfected tumor cells, U , grow at rate r. Virions, V , infect uninfected cells at a frequency-dependent rate with infection parameter β. T cells, T , at the tumor site kill uninfected cells at a frequency-dependent rate with coefficient k(I), which is a function of the number of infected cells, I.
• In (2), βU V /N is the rate that tumor cells become newly infected. Infected cells lyse at rate δ I and, like uninfected cells, are killed by T cells at a frequency-dependent rate. • In (3), u V (t) is the rate that new virions are injected into the system at time t, which is determined by the treatment protocol. Upon lysis of each infected cell, α new virions are released and free virions are cleared at rate δ V . • In (4), 4-1BBL provides a co-stimulatory signal to activate T cells at rate s T (I). Additional T cells are also activated by APCs in the presence of IL-12 and by injected DCs at rates χ A and χ D , respectively. T cells die at rate δ T . • In (5), IL-12 stimulates recruitment of mature, tumor-antigen-presenting APCs at rate s A (I). APCs then migrate to tumor-draining lymph nodes where they stimulate T cells that then migrate back to the tumor. Tumor-antigenpresenting APCs die or turnover at rate δ A . • In (6), u D (t) is the rate that tumor-antigen-presenting DCs are injected into the system at time t, based on the treatment protocol. DCs migrate to the lymph nodes and stimulate T cells in a similar manner to endogenous APCs. Tumor-antigen-presenting DCs die or turnover at rate δ D .
The major assumptions and conditions of the model are as follows: 1. T cell recruitment is promoted by 4-1BBL, which is produced by virus-infected tumor cells. We assume that the increase in T cell recruitment is a function of the amount of 4-1BBL at the tumor site, which we assume is proportional to the number of infected cells, so the T cell supply rate, s T (I), is a function of I. For simplicity, we suppose s T (I) = c T I where c T is constant. 2. APC recruitment is promoted by IL-12, which is produced by virus-infected tumor cells. As above, we assume that the increase in APC recruitment is a function of the amount of IL-12 at the tumor site, which we assume is proportional to the number of infected cells, so the APC supply rate, s A (I), is a function of I, and s A (I) = c A I. 3. T cell cytotoxicity, i.e. the T cell killing rate k(I), is influenced by both cytokines, 4-1BBL and IL-12. We assume the simple form k(I) = k 0 + c k I. The constant k 0 is the base T cell killing rate that can be increased linearly by the presence of cytokine-producing infected cells. 4. The variables T and A represent T cells and APCs recruited by virus-mediated 4-1BBL or IL-12. We implicitly assume that the baseline immune response is already factored into the net tumor growth rates. 5. While IL-12 can stimulate the innate immune response (e.g. natural killer cells), we only consider its effect on the cytotoxic T cell response [11], and we do not explicitly model an immune response, e.g. antibody response, to the virus itself. 6. In [18], tumor size is measured by volume (mm 3 ). In our model, it is measured by number of cells. To convert between volume and cell number, we assume tumor cells have a diameter of approximately 0.01 mm [2,6,29], and hence occupy a volume of order 10 −6 mm 3 , meaning that 1 mm 3 contains on the order of 10 6 tumor cells. The same system of equations could be obtained without such scaling, but the scaling is done to more easily describe our results in units consistent with those of the experiments. 7. We assume that all populations start at size 0, except for the initial population of uninfected tumor cells, U (0), which we fit to the experimental data. 2.2. Model fitting. In contrast to our earlier work [22], where we only considered treatment using immunostimulatory adenoviruses, in this paper, we used the extended model (1)-(6) to simultaneously fit the data of Figures 2a and 4a of Huang et al. [18], which corresponds to treatment using immunostimulatory OV alone and combining OV and DC vaccines.
2.2.1. Estimated parameters. Before fitting the model to data, we estimated seven of the parameters based on other information from experimental literature. We estimated that oncolytic adenovirus produce a viral burst size of α = 3000, which falls within the range measured by [7]. The time for infected cells to undergo lysis is on the order of 1 day, δ I = 1/day, which corresponds to an average lysis time of one day [4,16,17,19]. Ninety percent of viruses decay or exit the tumor site in one day [44]. Accordingly, we estimated a viral decay rate of δ V = − log(10%) = 2.3/day, which also agrees with the decay rates of Li et al. and Wang et al. [28,40]. We estimated that T cells have a half-life of 48 hours, so that δ T = 0.35/day as shown in [8].
In the case of T cell killing of tumor cells in mice, de Pillis et al. estimated a maximum fractional kill rate of 1.43/day [10]. We note that all experiments considered in this paper were performed on mice, so we do not consider the estimates obtained from human patients. Although their expression for the fractional kill rate in [10] differs from our expression in (1) and (2), we use their parameter as a ballpark value. Using this information, we estimated a base T cell kill rate of k 0 = 2/day. In fact, we tested the model over a range of values and found that base kill rates varying from 1.5 to 2.5 all resulted in reasonable fits of the model to data.
We assumed that APCs die or exit the system at a similar rate as T cells, yielding an APC death rate of δ A = 0.35/day. We used the same estimate δ D = 0.35/day for the decay rate of injected DCs. We also assumed that APCs activate a T cell in an average time of 1 day, so that χ A = 1/day [38,39]. These parameter estimates are summarized in Table 1.

Fit parameters.
To estimate the rest of the model parameters, we fit numerical solutions of the model (1) Tumor growth was recorded daily and represents the average of 7-9 murine models. All data sets contain measurements of murine B16-F10 melanoma tumor volumes with respect to time [18]. For the control case, we used the combined data from Figures 2a and 4a from [18], a data set representing 15 murine tumors from day 0 to 10 and 9 murine tumors on days 11 and 12. We fit the simplest possible model, exponential growth, for tumor growth in the absence of treatment. All parameters apart from the net tumor growth rate, r, were set to 0 (see the PBS column in Table 1). We estimated r using a weighted linear least-squares fit of the log of the tumor data. For our weighting, we normalized distances by the standard error at each time point. We obtained the fit r = 0.32/day and retained this estimate for the remaining models. See Figure 2.
For the remaining data sets, model parameters were fit sequentially and hierarchically by estimating parameters for submodels and retaining these values as the models became increasingly more advanced. We numerically simulated the ODE system (1)-(6) using the function 'ode45' in MATLAB 2014a (The MathWorks, Inc.; Natick, MA, USA) and estimated the remaining parameters by obtaining the weighted-least-squares fit of the model solution to data using the Levenberg-Marquadt algorithm, provided by the MATLAB function 'lsqnonlin.' As before, we weighted distances by the reciprocal of the standard error at each time point.
Unlike in the control case, we used the rectangular least-squares fit rather than the log-least squares fit. The primary reason for this change was that the high variance among individual mice under treatment (vii) using combination adenovirus (iv) Data and model simulation for treatment with Ad/IL-12 using doses of 5 × 10 9 virions on days 0, 2, and 4. Model parameters for cases (ii), (iii), and (iv) were estimated using weighted least squares fitting. Virus doses were taken from experimental protocols of [18].
and DCs caused the geometric mean of tumor sizes to deviate well below the arithmetic mean. In addition, it proved harder to fit the geometric mean in this case, so for consistency we chose to use arithmetic means and rectangular least-squares fitting for all treatment cases. Apart from case (vii), this choice made very little difference in fitting the other treatment cases (i)-(vi).
To fit the data corresponding to treatment (ii) with adenovirus alone, we used a model that excludes the infection-induced immune response. This submodel includes the additional parameters α (viral burst size), δ I (lysis rate of infected cells), and δ V (viral decay rate), as estimated from the literature (see Table 1). We fit the infection parameter β to match the data for adenovirus without cytokines and obtained the estimate β = 0.001/day. For the data sets corresponding to treatments (iii) and (iv) with Ad/4-1BBL and Ad/IL-12, we assumed that the viruses recruit T cells or APCs, respectively, and both increase the T cell killing rate. Thus, for both viruses, we fit the T cell killing rate constant, c k . For Ad/4-1BBL, we additionally fit the T cell recruitment parameter, c T , and for Ad/IL-12, we additionally fit the APC recruitment parameter, c A . The parameter estimates from our fits are shown in Table 1 and plots of experimental data and model simulations are shown in Figure 2.
For treatment (v), we simultaneously fit the model to the two data sets for Ad/4-1BBL/IL-12. In data set 1, doses of 5 × 10 9 virions were injected on days 0, 2, and 4 [18, Fig. 2a], and in data set 2, doses of 2.5 × 10 9 virions were injected on days 0, 2, and 4 [18, Fig. 4a]. Parameter estimates from these fits are shown in Table 2 and plots of experimental data and model simulations are shown in  We also simultaneously fit data for treatments (vi) and (vii). In this case, we inherited all parameters from the previous fit for treatment (v) and only fit the parameter χ D , the rate of T cell stimulation by DCs from injected DC vaccines. For treatment (vi), 10 6 DCs were injected on days 1, 3, and 5. For treatment (vii), 2.5 × 10 9 particles of Ad/4-1BBL/IL-12 virus were injected on days 0, 2, and 4, and 10 6 DCs were injected on days 1, 3, and 5. Parameter estimates from these fits are shown in Table 2     Base T cell killing rate (day −1 ) - T cell stimulation rate by APCs (day −1 ) ---1  Table 2. Parameter estimates used in the model equations (1)- (6). Parameters in the table were estimated from data using weighted least-squares fitting. Other parameters were fixed at the values shown in the bottom half of Table 1.

Results.
3.1. The fit. Visual inspection suggests that the nonlinear fitting methods used to produce parameters worked well with simpler models (Figure 2). Increasingly more complex models are progressively less accurate, as one would expect, but still reasonably well fit. As we can see from Figure 4(vii), the model simulation for combination treatment (vii) using viruses and DC vaccines does not appear to fit the data as closely as the other cases. However, if we view the experimentally measured and simulated tumor plots on a logarithmic scale, we see that the model still closely captures much of the behavior of the experimental data. Figure 5 shows how the model follows the rising, falling, and rising behavior of the data and captures the long-term exponential rate of increase of the tumor volume after day 15. Pearson's r coefficient was calculated to give a more rigorous measurement of fit. Agreement between the model and the data is found when Pearson's r is near 1. For all models, Pearson's r was over 0.9 (Table 3), and over 0.98 for all but the most complicated model (which had r = 0.92), suggesting good fits for all of the models.  Huang et al. used a DC dose of 1 × 10 6 DCs. Using this dose, tumor elimination did not occur. To determine if there was a bifurcation value (of DC dose), that when crossed, would cause the tumor-free state to be approached, we simulated treatment varying DC doses. Figure 6 shows end tumor volume versus DC dose for doses of DC ranging between 0.94 × 10 6 DCs and 1.04 × 10 6 DCs. Tumor volume is below the threshold (1 × 10 −6 mm 3 ) for DC doses greater than 1.019 × 10 6 DCs. The DC dose chosen by Huang et al. was near the bifurcation value found here. In fact, in the experimental data of Huang et al., three out of eight mice did achieve complete remission, lending support to our result that the chosen DC dose was near a bifurcation value [18].

Alternating OV and DC injections.
Huang et al. sought to discover if combination therapy with OV and DC injections was more effective than either treatment alone. The treatment strategy investigated was three alternations of OV injections with DC injections (6 total injections). They found that the combination therapy was more effective than three doses of OV given at days 0, 2 and 4, and also more effective than three doses of DC injected at days 1, 3, and 5.
To investigate whether alternating treatments was most effective, we simulated all twenty combinations of three OV and three DC injections. One strategy (OV-OV-OV-DC-DC-DC) caused elimination of the tumor (volume < 10 −6 mm 3 ). For the other nineteen strategies, the tumor was not eliminated and tumor burden remained between 60 mm 3 and 1,100 mm 3 at 30 days. There was no apparent order to which strategy did best or worst ( Table 4 Table 4. Treatment order at days 0, 1, 2, 3, 4, 5 and resulting tumor volume at 30 days. OV doses were 2.5 × 10 9 virions and DC doses were 1 × 10 6 DCs. To further investigate optimal treatment strategies, we simulated all sixty-four possible six-injection strategies. Strategies containing more than three OV injections did the worst (largest end tumor volume), but not in any apparent order. For example, the worst ten strategies contained two DC injections while the 22nd worst was all OV injections.
In sharp contrast, those strategies containing more than three DC injections eliminated the tumor altogether (defined as volume < 1 × 10 −6 mm 3 ). Given that one dose of DC results in a larger end tumor volume than one dose of OV (Figure 7), we found this result surprising. With one injection, OV may have a stronger effect than DC injections because of the combined effects from IL-12 and 4-1BBL. However, injected DC cells have a larger effect on T cell recruitment than APCs in the presence of IL-12 (χ D > χ A ), which may explain why multiple DC injections cause a larger reduction in tumor burden.
Additionally, treating with six DC injections caused tumor remission suggesting that a strong reduction in tumor burden could be achieved by increasing DC treatment even in the absence of OV. In contrast, when the dose of OV alone was doubled on days 0, 2 and 4 from 2.5 × 10 9 to 5 × 10 9 virions, tumor burden was reduced but still not as effectively as with combination therapy of OV and DC injections [18,Figs. 2a,4a]. Therefore, DC injections enhance the effectiveness of OV, but while DCs exhibit strong anti-tumor activity without OVs, the complementary and synergistic role of OV is a useful clinical observation.

3.2.3.
Varying dose size. In the previous section, the results showed that strategies containing more than three DC injections performed better (less tumor burden after 30 days) than those containing three or fewer DC injections. To determine the effect of the size of the injections on this result, we simulated treatment with OV and DC, varying the injections amounts of each.
Using the strategy of Huang et al., six injections of OV and DC were simulated in the order OV-DC-OV-DC-OV-DC at times 0, 1, 2, 3, 4, and 5 days respectively. The goal in varying the doses was to find the smallest doses of OV and DC that would eliminate the tumor after 30 days (defined by tumor volume less than 10 −6 mm 3 ) and then to see how the doses differed based on the order of injections. To determine the range of the doses, we first simulated treatment with OV alone at days 0, 2, and 4 days and found the smallest dose that would eliminate the tumor in the absence of DC injections. The minimum dose for OV injections alone was 6.3692 × 10 9 virions. Similarly, we simulated treatment with DC injections at days 1, 3, and 5, in the absence of OV injections, looking for the minimum dose necessary to eliminate the tumor at day 30 and found the dose to be 2.066 × 10 6 DCs. Next OV doses were varied from 0 to 6.3 × 10 9 virions and the smallest DC dose needed to eliminate the tumor was found. These bifurcation values are shown in Figure 8.
Next, the same simulations were performed for different combinations of injections. The results are shown in Figure 8. Interestingly, the curves intersect near the  suggesting that strategies using these doses are similar. Note also that as the doses change, the ranked efficacy of the strategies completely reverse, showing that the goodness of the strategy is dependent on the dose of OV and DC. The strategies show the biggest differences when one dose is small suggesting that the ordering of the injections matter. This is notable when OV or DC doses are near zero, so that the strategies differ only by order of injection.
When DC injections are lumped together and given first, OV has a smaller effect with small DC doses requiring a huge increase of OV doses to eliminate the tumor.
In contrast, when OV injections are lumped together and given first, DC injections can be altered between 0 and 3.6 × 10 6 cells to cause tumor elimination.

3.2.4.
Varying time between injections and number of doses. We also investigated the effect of varying the time between consecutive doses and the number of doses. For simplicity, we only considered cases where all doses are evenly-spaced in time.
Based on the strategy of Huang et al., we considered an alternative pattern of OV followed by DC injections [18]. We considered sequences of 2, 4, 6, 8, and 10 doses spaced at time intervals between 0.01 and 3 days. The OV dose was set to 2.5 × 10 9 virions, and the DC dose was set to 1 × 10 6 DCs as in Huang et al. [18]. Figure 9 shows tumor volumes at day 30 versus the time between injections for various numbers of doses. Time between injections (day) Figure 9. Plots of tumor volumes at day 30 versus times between injections for treatments of n = 2, 4, 6, 8, and 10 doses of alternating OV and DCs. OV dose was set to 2.5 × 10 9 virions, and DC dose was set to 1 × 10 6 DCs as in Huang et al. [18].
As expected, a greater number of doses always results in greater reduction of tumor volume, so the curve corresponding to 10 doses is below the blue curve corresponding to 8 doses and so on. We also observed that a stronger response to doses given closer together, converging to a heavy dose all at once. In fact, a heavy dose given all at once even leads to a lower tumor burden on day 30 than a drawn out treatment strategy over several days. This result might occur because a heavy dose given all at once or almost all at once leads to a heavier initial virus infection and a stronger initial immune response. As with all treatments, the maximum intensity and number of doses that can be used will be limited by toxicity considerations, which we will consider in a future work [20,21].
We compare the results in Figure 9 to the following two alternative treatment combinations: (1) alternating DC and OV, and (2) a sequence of OV doses followed by an equal number of DC doses. The results of the simulations are shown in Figure 10.
As we can see from Figure 10 OV dose was set to 2.5 × 10 9 virions, and DC dose was set to 1 × 10 6 DCs as in Huang et al. [18].
in this case the two patterns are effectively the same. On the other hand, from Figure 10(right), for a fixed number of doses, a sequence of OV doses followed by an equal number of DC doses always does better than the pattern of alternating OV and DC doses, except when the patterns are the same for n = 2. From these results, it seems that the ranking implied by Table 4, that a sequence of OV followed by a sequence of DC consistently performs better than the original strategy of alternating OV and DC, remains the same even when varying the numbers of doses and the time intervals between them.
These findings strongly suggest that it is more effective to treat a tumor with immunostimulatory oncolytic viruses first to induce a potent and expanding infection without enhancing the immune response, which would reduce the infection, and then follow up with a sequence of DCs to stimulate a subsequent anti-tumor immune response that will temporarily work synergistically with the viral infection. The idea of temporally separating viral treatment and DC vaccines is a different strategy than what was tried in the experiments of Huang et al. [18] and would be an interesting alternative to consider in future in vivo investigations. 4. Summary. Experimental models and recent clinical trials have shown oncolytic virotherapy to be a viable, safe treatment for cancer. The inability of virotherapy to completely eliminate tumors seems to be the current major limitation, giving rise to the investigation of combination treatment protocols. In a recent work, Huang et al. showed that combination therapy consisting of alternating injections of immunostimulatory oncolytic virus (OV) with dendritic cell (DC) vaccines was more effective than either treatment in isolation [18].
Alternating OV injections with DC injections at the particular dose that Huang et al. used is one of many possible treatment strategies for combining OV and DC therapies. Experimentally exploring other strategies is costly, requiring both time and physical resources, and even well-funded studies can explore only a small range of possibilities. In this paper, alternative strategies were explored using a mathematical model, a scientific tool that allows for an easier, less-costly, and more thorough investigation of treatment regimens.
We developed a mathematical model consisting of ordinary differential equations to model tumor growth based on the murine models of Huang et al. [18]. We chose our model to be simple enough to estimate realistic parameter values and to be able to hierarchically fit the model based on the data of Huang et al., that utilized increasingly more complicated treatments. While it is possible that the hierarchical fitting may result in sub-optimal parameters, Pearson's r coefficients showed that our model fit the data well, with all coefficients being over 0.92. We then used this model to make predictions about other treatment strategies. Although the model well-describes the data for all treatment protocols, it is still important to point out that the conclusions we draw from our analysis may indeed be model-dependent. Therefore, in the future, we will test how sensitive the predictions are to the design of the model, as well as test sensitivity to various model parameters.
We made the following observations: 1) a larger dose of DC could cause complete elimination of the tumor; 2) combinations other than alternating OV-DC, the approach chosen by Huang et al. would further reduce tumor burden or eliminate the tumor completely, with a protocol of all OVs followed by all DCs being most effective; 3) the efficacy of a particular protocol (ordering of OVs and DCs) depends on the dose of each, and different strategies would be more effective based on which doses were chosen; 4) the efficacy of treatment would increase as doses were given closer together; and 5) the efficacy of a dosing strategy depends on the temporal spacing of the doses. Again, while some of the particular details may be modeldependent, the strong model fit to the data lends credence to these conclusions.
Addressing the first hypothesis, results from simulations showed that the DC dose chosen by Huang et al. was near the bifurcation point such that when increased from 1 × 10 6 DCs to 1.02 × 10 6 DCs, trajectories approached the tumor-free state. This result was particular to the treatment strategy chosen by Huang et al. (alternating OV and DC injections) and to the model developed. This result suggests that the dose of DCs can be increased to cause complete elimination of the tumor (volume of less than 10 −6 at 30 days). However, we note that increasing doses can lead to systemic toxicity.
Secondly, altering OV and DC injections doses, as done by Huang et al., does not seem to be the most effective treatment strategy. We found that this strategy ranked 15th among the twenty possible strategies containing three OV doses and three DC doses for doses of 2.5×10 9 virions and 1×10 6 DCs. Moreover, simulations support the notion that separating OV and DC doses is more effective for reducing tumor burden than alternating these injections. Furthermore, whether one should inject OV or DC first is dependent on the chosen doses of each. In general, the injection sequence can be chosen chosen such that adverse effects are minimized.
Upon further investigation, the choice of strategy depends on the dose size of OV and DC injections. Interestingly, all strategies eliminated the tumor for doses just slightly larger than those given by Huang et al. [18]. As either DC or OV dose was increased, a smaller dose of the other was needed to cause tumor elimination, but the amount of the increase or decrease varied with the chosen strategy.
Finally, results from simulations show that injections of OVs and DCs are most effective when given closer together, and in fact the strongest response occurs when all injections are given at the beginning of the trial. Giving a larger number of doses increases the efficacy of treatment, but there appears to be an asymptotic level of efficacy approached for any given treatment as the dose is increased. As mentioned in Section 3.2.4, the optimal dose magnitude and number will by restricted by considering the expected toxicity of the treatment strategy, an extension that we intend to consider in a future work [20,21].
Overall, results from our simulations show that tumor elimination can be achieved by finely tuning doses of OV and DC injections and by choosing strategies based on the particular dose chosen. Oncolytic virotherapy in combination with other treatment protocols is progressing in efficacy. Our work shows how the chosen treatment strategy and doses given can impact efficacy of treatment on tumor reduction and elimination.