Friday Greetings,
I'm working on an ODE system as a basic model for cancer treatment:
$\begin{array}{c} \frac{{dN}}{{dt}} = kN\left( {1 - \frac{N}{\alpha }} \right) - \mu CN\\ \frac{{dC}}{{dt}} = d\left( t \right) - \lambda C - \gamma CN\\ N\left( {{t_0}} \right) = {N_0}\\ C\left( {{t_0}} \right) = {C_0} \end{array}$
where
$d\left( t \right) = \left\{ {\begin{array}{*{20}{c}} {{\rm{dosage}}}&{}&{n\left( {T - 1} \right) < t < n\left( {T - \tau } \right)}\\ 0&{}&{n\left( {T - \tau } \right) < t < nT} \end{array}} \right.\,\,\,\,{\rm{with}}\,\,\,\,n \in \left\{ {1,\,2,\,3,\,...} \right\}$
Ultimately I want to produce an accurate graph of the treatment and response. Here is my code:
k = 1; alpha = 1; mu = 1; lambda = 1; gamma = 0.5; tnot = 0;
nnot = 0.1; cnot = 0; tmax = 84;
nmin = -0.1; nmax = 1.15; cmin = -0.1; cmax = 1.15;
doseinterval = 10/24; wholeinterval = 2; numberofdoses = 30;
applieddosingPeriodic[t_] := If[t < wholeinterval*numberofdoses, If[Mod[t, wholeinterval] < doseinterval, dosage, 0], 0];
paramSolnCancerTreatPeriodic = ParametricNDSolveValue[{n'[t] == k*n[t]*(1 - n[t]/alpha) - mu*c[t]*n[t], c'[t] == applieddosingPeriodic[t] - lambda*c[t] - gamma*c[t]*n[t], n[tnot] == nnot, c[tnot] == cnot}, {n, c}, {t, tnot, tmax}, {dosage}, MaxStepSize -> 0.01];
nPlotCancerTreatPeriodic[dosage_] := Plot[paramSolnCancerTreatPeriodic[dosage][[1]][t], {t, tnot, tmax}, PlotStyle -> {Thick, Blue}, AxesLabel -> {Style[t, 16], Row[{Style["N", 16, Blue], Style[", ", 16, Black], Style["C", 16, Red]}]}, PlotLabel -> Style["Cancer Treatment with Periodic Dosing", 16], PlotRange -> {{tnot - 2, tmax + 2}, {nmin, nmax}}, Epilog -> {Black, PointSize[Large], Point[{{tnot, nnot}, {tnot, cnot}, {tmax, paramSolnCancerTreatPeriodic[dosage][[1]][tmax]}, {tmax, paramSolnCancerTreatPeriodic[dosage][[2]][tmax]}}], Green, PointSize[Medium], Point[{tnot, nnot}], Point[{tnot, cnot}], Red, Point[{tmax, paramSolnCancerTreatPeriodic[dosage][[1]][tmax]}], Point[{tmax, paramSolnCancerTreatPeriodic[dosage][[2]][tmax]}]}];
cPlotCancerTreatPeriodic[dosage_] := Plot[paramSolnCancerTreatPeriodic[dosage][[2]][t], {t, tnot, tmax}, PlotStyle -> {Thick, Red}, AxesLabel -> {Style[t, 16], Row[{Style["N", 16, Blue], Style[", ", 16, Black], Style["C", 16, Red]}]}, PlotLabel -> Style["Cancer Treatment with Periodic Dosing", 16], PlotRange -> {{tnot - 2, tmax + 2}, {cmin, cmax}}, Epilog -> {Black, PointSize[Large], Point[{{tnot, nnot}, {tnot, cnot}, {tmax, paramSolnCancerTreatPeriodic[dosage][[1]][tmax]}, {tmax, paramSolnCancerTreatPeriodic[dosage][[2]][tmax]}}], Green, PointSize[Medium], Point[{tnot, nnot}], Point[{tnot, cnot}], Red, Point[{tmax, paramSolnCancerTreatPeriodic[dosage][[1]][tmax]}], Point[{tmax, paramSolnCancerTreatPeriodic[dosage][[2]][tmax]}]}];
frameCancerTreatPeriodic[dosage_] := Show[nPlotCancerTreatPeriodic[dosage], cPlotCancerTreatPeriodic[dosage], ImageSize -> 1040];
frameCancerTreatPeriodic[1.0]
The last command produces the following image:
graphical output with aberration http://www.biomathdynamics.com/mmartin/cancerPeriodicTreatError.jpg
Notice the slight aberration in the blue response. The adaptive stepsize in ParametricNDSolveValue (or, ultimately, NDSolve) is too large for one part of the graph (around a time of 25 to 30). I've tried tweaking several of the different options/parameters for NDSolve, but nothing seems to efficiently work. I have the sense that it's something simple that I'm just missing. I would appreciate your input and guidance on that as I clean up my modeling efforts in this domain. Ultimately, I am making an animated gif as the dosage parameter is varied and, depending on the dosage, the aberration may be evident or not.
Thank you in advance for any help/advise.
Regards, Mike
PlotPoints
? (Welcome to Mathematica.SE!) $\endgroup$PlotPoints
isn't a listed parameter for eitherNDSolve
orParametricNDSolveValue
. I'm using v9 but have v10 on the way -- maybe it's included there. Your suggestion makes sense, though. $\endgroup$Plot
, notNDSolve
. $\endgroup$