﻿ Optimal Control of Cancer Growth

Applied Mathematics
Vol.10 No.04(2019), Article ID:91699,23 pages
10.4236/am.2019.104014

Optimal Control of Cancer Growth

Jens Christian Larsen

Vanløse Alle 50 2. mf. tv, 2720 Vanløse, Copenhagen, Denmark    Received: February 15, 2019; Accepted: April 8, 2019; Published: April 11, 2019

ABSTRACT

The purpose of the present paper is to apply the Pontryagin Minimum Principle to mathematical models of cancer growth. In  , I presented a discrete affine model T of cancer growth in the variables C for cancer, GF for growth factors and GI for growth inhibitors. One can sometimes find an affine vector field X on ${ℝ}^{3}$ whose time one map is T. It is to this vector field we apply the Pontryagin Minimum Principle. We also apply the Discrete Pontryagin Minimum Principle to the model T. So we prove that maximal chemo therapy can be optimal and also that it might not depending on the spectral properties of the matrix A, (see below). In section five we determine an optimal strategy for chemo or immune therapy.

Keywords:

Cancer, Models of Cancer Growth, Pontryagin Minimum Principle 1. Introduction

Our reference to optimal control theory of ODEs is  and to control theory of discrete systems . There is a review of papers in optimal control theory applied to cancer . The model we consider here is from  and is defined by

$A=\left(\begin{array}{ccc}1+\gamma & \alpha & \beta \\ \delta & 1+{\mu }_{F}& 0\\ \sigma & 0& 1+{\mu }_{I}\end{array}\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}T\left(y\right)=Ay+g$ (1)

where $y={\left(C,GF,GI\right)}^{\text{T}}$,$g={\left({g}_{C},{g}_{F},{g}_{I}\right)}^{\text{T}}\in {ℝ}^{3}$, and T denotes transpose. Here $\alpha \in {ℝ}_{+}$,$\beta \in {ℝ}_{-}$,$\delta ,\sigma \in ℝ$,${\mu }_{F},{\mu }_{I}<0$. We assume that $\alpha \delta +\beta \sigma \ne 0$.

The matrix A has characteristic polynomial $-p\left(\lambda \right)$ where

$p\left(\lambda \right)={\lambda }^{3}-{\lambda }^{2}\left(3+\gamma +{\mu }_{F}+{\mu }_{I}\right)-\lambda \left(\alpha \delta +\beta \sigma -\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)$ (2)

$\text{ }-\left(1+\gamma \right)\left(2+{\mu }_{F}+{\mu }_{I}\right)+\alpha \delta \left(1+{\mu }_{I}\right)+\beta \sigma \left(1+{\mu }_{F}\right)-\left(1+\gamma \right)\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)$ (3)

and this polynomial can have (i) three real roots or (ii) one real root and two imaginary roots. It turns out that asking $\mu ={\mu }_{F}={\mu }_{I}$ simplifies matters considerably. Then

$p\left(\lambda \right)=-\left(1+\mu -\lambda \right)\left({\lambda }^{2}-\left(2+\gamma +\mu \right)\lambda +\left(1+\gamma \right)\left(1+\mu \right)-\alpha \delta -\beta \sigma \right)$ (4)

So the eigenvalues are $1+\mu$ and

${\lambda }_{±}=1+\frac{\gamma +\mu }{2}±\frac{1}{2}\sqrt{{\left(\gamma -\mu \right)}^{2}+4\left(\alpha \delta +\beta \sigma \right)}$ (5)

We assume that $1+\mu >0,\mu <0$. In case (i) if the eigenvalues of A, ${\lambda }_{+},{\lambda }_{-},\stackrel{˜}{\lambda }$ are positive and distinct and we assume this, then you can find an affine vector field X on ${ℝ}^{3}$ such that the time one map is d

${\Phi }_{1}^{X}=T$ (6)

see  and below. In case (ii) if the eigenvalues of A are $1+\mu ,a±ib,a>0,1+\mu >0,b\ne 0$ and we assume this, you can also find an affine vector field X on ${ℝ}^{3}$ such that

${\Phi }_{1}^{X}=T$ (7)

In case (i) define a matrix of eigenvectors

$D=\left(\begin{array}{ccc}1+\mu -{\lambda }_{+}& 1+\mu -{\lambda }_{-}& 0\\ -\delta & -\delta & -\beta \\ -\sigma & -\sigma & \alpha \end{array}\right)$ (8)

and in case (ii) define

$U=\left(\begin{array}{ccc}1+\mu -a& -b& 0\\ -\delta & 0& -\beta \\ -\sigma & 0& \alpha \end{array}\right)$ (9)

Then in case (i)

$\Lambda ={D}^{-1}AD=\left(\begin{array}{ccc}{\lambda }_{+}& 0& 0\\ 0& {\lambda }_{-}& 0\\ 0& 0& \stackrel{˜}{\lambda }\end{array}\right)$ (10)

if the eigenvalues are distinct and positive and in case (ii)

$\Lambda ={U}^{-1}AU=\left(\begin{array}{ccc}a& b& 0\\ -b& a& 0\\ 0& 0& 1+\mu \end{array}\right)$ (11)

see . Here ${D}^{-1}$ denotes the inverse to D. To find an X in case (i) define when the eigenvalues are real, positive and distinct, an affine vector field

$Y\left(x\right)=\left(\begin{array}{ccc}\mathrm{ln}{\lambda }_{+}& 0& 0\\ 0& \mathrm{ln}{\lambda }_{-}& 0\\ 0& 0& \mathrm{ln}\stackrel{˜}{\lambda }\end{array}\right)x+d$ (12)

Then with the right choice of d

${\Phi }_{1}^{Y}\left(x\right)={D}^{-1}ADx+{D}^{-1}g$ (13)

hence if we let

$X\left(y\right)=DY\left({D}^{-1}\left(y\right)\right)$ (14)

we get

${\Phi }_{1}^{X}=T$ (15)

And to find an X in case (ii) when $a>0,b\ne 0,1+\mu >0$, define

$Y\left(x\right)=\left(\begin{array}{ccc}{a}_{1}& {b}_{1}& 0\\ -{b}_{1}& {a}_{1}& 0\\ 0& 0& \mathrm{ln}\left(1+\mu \right)\end{array}\right)+\stackrel{˜}{d}$ (16)

Then with the right choice of ${a}_{1},{b}_{1}\in ℝ,\stackrel{˜}{d}\in {ℝ}^{3}$ we get

${\Phi }_{1}^{Y}\left(x\right)={U}^{-1}AU\left(x\right)+{U}^{-1}g$ (17)

hence with

$X\left(y\right)=UY\left({U}^{-1}\left(y\right)\right)$ (18)

we have that the time one map is

${\Phi }_{1}^{X}\left(y\right)=T\left(y\right)$ (19)

See  for details of the above and also below.

In section two we solve the problem:

minmize $C\left(T\right)$,$T>0$ (fixed) subject to

${\left(\begin{array}{c}C\\ GF\\ GI\end{array}\right)}^{\prime }=X\left(\begin{array}{c}C\\ GF\\ GI\end{array}\right)+{e}_{3}u\left(t\right)$ (20)

$u\left(t\right)\in \left[0,{g}_{I}^{0}\right],{g}_{I}^{0}>0$, by first solving it for Y and then infer the solution for X. Here

${e}_{3}=B={\left(0,0,1\right)}^{\text{T}}$ (21)

and $u\left(t\right)$ is piecewise continuous. In section three we apply the discrete Pontryagin minimum principle to the difference equation

${x}_{k+1}=A{x}_{k}+B{u}_{k}+g$

where

${x}_{k}={\left(C\left(k\right),GF\left(k\right),GI\left(k\right)\right)}^{\text{T}}$ (22)

${u}_{k}\in \left[0,{g}_{I}^{0}\right],{g}_{I}^{0}>0$. $k=0,1,\cdots ,N-1$, where $N\in ℕ,N\ge 2,{x}_{0}=x$, with the objective to minimize $C\left(N\right)$. There are again two cases to consider (i) and (ii) above. If $\mu ={\mu }_{F}={\mu }_{I}$ and the eigenvalues are positive and distinct, maximal chemo therapy is optimal. But in case (ii) it is not always optimal.

When ${\mu }_{F}\ne {\mu }_{I}$ and the eigenvalues are real and distinct, we produce a counter example to maximal chemo therapy being optimal, see section four. Some solid tumors grow like Gompertz functions, see . There are several important monographs in mathematics and medicine, see  - .  -  are my latest papers on cancer and mathematics.  is our reference to roots of cubic polynomials.  proves continuous dependence of roots of a polynomial on the coefficients of the polynomial.

In section five we consider optimality of the discrete model T when ${\mu }_{F}\ne {\mu }_{I}$. Here we also determine optimal control of the map T.

2. Optimal Control of X

The purpose of this section is to minimize $C\left(T\right)$, subject to

${\left(\begin{array}{c}C\\ GF\\ GI\end{array}\right)}^{\prime }=X\left(\begin{array}{c}C\\ GF\\ GI\end{array}\right)+{e}_{3}u\left(t\right)$ (23)

$T>0$ fixed and with $\mu ={\mu }_{F}={\mu }_{I}$. Let us consider (ii) first. We assume that there are a real eigenvalue $1+\mu$ and two imaginary eigenvalues $a±ib,a>0,1+\mu >0,b\ne 0$.

Now define the two by two matrix

${L}_{1}=\left(\begin{array}{cc}{a}_{1}& {b}_{1}\\ -{b}_{1}& {a}_{1}\end{array}\right)$ (24)

where

${b}_{1}={\mathrm{tan}}^{-1}\left(\frac{b}{a}\right)$ (25)

${a}_{1}=\mathrm{ln}\sqrt{{a}^{2}+{b}^{2}}$ (26)

This will imply, that

$\mathrm{exp}\left(\begin{array}{cc}{a}_{1}& {b}_{1}\\ -{b}_{1}& {a}_{1}\end{array}\right)=\left(\begin{array}{cc}a& b\\ -b& a\end{array}\right)=L$ (27)

Also let

$\stackrel{˜}{B}=\left(\begin{array}{ccc}{a}_{1}& {b}_{1}& 0\\ -{b}_{1}& {a}_{1}& 0\\ 0& 0& \mathrm{ln}\left(1+\mu \right)\end{array}\right)$ (28)

and define the vector field

$Y\left(x,v\right)=\stackrel{˜}{B}x+d+{U}^{-1}{e}_{3}v=\stackrel{˜}{Y}\left(x\right)+{U}^{-1}{e}_{3}v$ (29)

which is affine when $v=0$, where ${e}_{1},{e}_{2},{e}_{3}$ is the canonical basis in ${ℝ}^{3}$. Also $x,d\in {ℝ}^{3},v\in \left[0,{g}_{I}^{0}\right],{g}_{I}^{0}>0$. Put

$\stackrel{˜}{X}\left(y\right)=U\stackrel{˜}{Y}{U}^{-1}\left(y\right)$ (30)

Let

$X\left(y,v\right)=\stackrel{˜}{X}\left(y\right)+{e}_{3}v=UY\left({U}^{-1}\left(y\right),v\right)$ (31)

Define ${d}_{1},{d}_{2},{d}_{3}$ by

${U}^{-1}{\left(g\right)}_{1,2}={L}_{1}^{-1}\left(L-\text{id}\right){d}_{1,2}$ (32)

and

$\frac{\mu }{\mathrm{ln}\left(1+\mu \right)}{d}_{3}={\left({U}^{-1}g\right)}_{3}$ (33)

Then

${\Phi }_{1}^{X}=T$ (34)

when $v=0$. To this vector field with $v\in \left[0,{g}_{I}^{0}\right]$ associate the Hamiltonian

$H\left(x,p,v,t\right)={p}^{\text{T}}Y\left(x,v\right)$ (35)

Then we have the adjoint equations

${{p}^{\prime }}_{1}=-\frac{\partial H}{\partial {x}_{1}}=-{p}_{1}{a}_{1}+{p}_{2}{b}_{1}$ (36)

${{p}^{\prime }}_{2}=-\frac{\partial H}{\partial {x}_{2}}=-{p}_{1}{b}_{1}-{p}_{2}{a}_{1}$ (37)

${{p}^{\prime }}_{3}=-\frac{\partial H}{\partial {x}_{3}}=-{p}_{3}\mathrm{ln}\left(1+\mu \right)$ (38)

So

${{p}^{\prime }}_{1,2}=\left(\begin{array}{cc}-{a}_{1}& {b}_{1}\\ -{b}_{1}& -{a}_{1}\end{array}\right)\left(\begin{array}{c}{p}_{1}\\ {p}_{2}\end{array}\right)$ (39)

which has flow

${p}_{1,2}\left(t\right)=\mathrm{exp}\left(-{a}_{1}t\right)\left(\begin{array}{cc}\mathrm{cos}\left({b}_{1}t\right)& \mathrm{sin}\left({b}_{1}t\right)\\ -\mathrm{sin}\left({b}_{1}t\right)& \mathrm{cos}\left({b}_{1}t\right)\end{array}\right){p}_{1,2}\left(0\right)$ (40)

Define

$d\left(t\right)={L}_{1}^{-1}\left(\mathrm{exp}\left({L}_{1}t\right)-\text{id}\right){d}_{1,2}$ (41)

Then the flow of Y is for $v=0$

${\Phi }^{Y}{\left(t,x\right)}_{1,2}=\mathrm{exp}\left({L}_{1}t\right)\left(\begin{array}{c}{x}_{1}\\ {x}_{2}\end{array}\right)+d\left(t\right)$ (42)

${\Phi }^{Y}{\left(t,x\right)}_{3}=\mathrm{exp}\left(\mathrm{ln}\left(1+\mu \right)t\right){x}_{3}+\frac{\mathrm{exp}\left(\mathrm{ln}\left(1+\mu \right)t\right)-1}{\mathrm{ln}\left(1+\mu \right)}{d}_{3}$ (43)

Define

${S}_{1}\left(x\right)=\left(1+\mu -a\right){x}_{1}-b{x}_{2}$ (44)

Then we have the transversality conditions

${p}_{1}\left(T\right)=\frac{\partial {S}_{1}}{\partial {x}_{1}}=1+\mu -a$ (45)

${p}_{2}\left(T\right)=\frac{\partial {S}_{1}}{\partial {x}_{2}}=-b$ (46)

${p}_{3}\left(T\right)=\frac{\partial {S}_{1}}{\partial {x}_{3}}=0$ (47)

This means that ${p}_{3}\left(t\right)=0$ and

${p}_{1,2}\left(T\right)=\mathrm{exp}\left(-{a}_{1}T\right)\left(\begin{array}{cc}\mathrm{cos}\left({b}_{1}T\right)& \mathrm{sin}\left({b}_{1}T\right)\\ -\mathrm{sin}\left({b}_{1}T\right)& \mathrm{cos}\left({b}_{1}T\right)\end{array}\right){p}_{1,2}\left(0\right)=\left(\begin{array}{c}1+\mu -a\\ -b\end{array}\right)$ (48)

The two by two matrix in this equation has inverse

$\left(\begin{array}{cc}\mathrm{cos}\left({b}_{1}T\right)& -\mathrm{sin}\left({b}_{1}T\right)\\ \mathrm{sin}\left({b}_{1}T\right)& \mathrm{cos}\left({b}_{1}T\right)\end{array}\right)$ (49)

So

${p}_{1,2}\left(0\right)=\mathrm{exp}\left({a}_{1}T\right)\left(\begin{array}{cc}\mathrm{cos}\left({b}_{1}T\right)& -\mathrm{sin}\left({b}_{1}T\right)\\ \mathrm{sin}\left({b}_{1}T\right)& \mathrm{cos}\left({b}_{1}T\right)\end{array}\right)\left(\begin{array}{c}1+\mu -a\\ -b\end{array}\right)$ (50)

But then we have

${p}_{1,2}\left(t\right)=\mathrm{exp}\left({a}_{1}\left(T-t\right)\right)\left(\begin{array}{cc}\mathrm{cos}\left({b}_{1}\left(t-T\right)\right)& \mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\\ -\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)& \mathrm{cos}\left({b}_{1}\left(t-T\right)\right)\end{array}\right)\left(\begin{array}{c}1+\mu -a\\ -b\end{array}\right)$ (51)

Notice that we have

${U}_{31}=b\beta$ (52)

${U}_{32}=\left(1+\mu -a\right)\beta$ (53)

If

$H\left({x}^{*}\left(t\right),p\left(t\right),{u}^{*}\left(t\right),t\right)\le H\left({x}^{*}\left(t\right),p\left(t\right),u,t\right)$ (54)

for all $u\in \mathbb{U}=\left[0,{g}_{I}^{0}\right]$, then $\left({x}^{*}\left(t\right),{u}^{*}\left(t\right)\right)$ is optimal, see Equations (55) to (60). But this amounts to the inequality

$\left({p}_{1}\left(t\right)b\beta {u}_{*}+{p}_{2}\left(t\right)\left(1+\mu -a\right)\beta {u}_{*}\right)\frac{1}{\mathrm{det}\left(U\right)}$ (55)

$=\mathrm{exp}\left({a}_{1}\left(T-t\right)\right)\left(\left(\mathrm{cos}\left({b}_{1}\left(t-T\right)\right)\left(1+\mu -a\right)-b\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\right)b\beta {u}_{*}$ (56)

$+\left(-\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\left(1+\mu -a\right)-b\mathrm{cos}\left({b}_{1}\left(t-T\right)\right)\right)\left(\left(1+\mu -a\right)\beta {u}_{*}\right)\right)\frac{1}{\mathrm{det}\left(U\right)}$ (57)

$=\mathrm{exp}\left({a}_{1}\left(T-t\right)\right)\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\left(-{b}^{2}-{\left(1+\mu -a\right)}^{2}\right)\beta {u}_{*}\frac{1}{\mathrm{det}\left(U\right)}$ (58)

$=\mathrm{exp}\left({a}_{1}\left(T-t\right)\right)\left(\alpha \delta +\beta \sigma \right)\beta {u}_{*}\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\frac{1}{\mathrm{det}\left(U\right)}$ (59)

$\le \mathrm{exp}\left({a}_{1}\left(T-t\right)\right)\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)\frac{\beta u}{-b}$ (60)

where $\text{det}\left(U\right)$ denotes the determinant of U. We have used that

$1+\mu -a=1+\mu -\left(1+\frac{\gamma +\mu }{2}\right)=\frac{\mu -\gamma }{2}$ (61)

and

${\lambda }_{+}=1+\frac{\gamma +\mu }{2}+\frac{1}{2}\sqrt{{\left(\gamma -\mu \right)}^{2}+4\left(\alpha \delta +\beta \sigma \right)}=a+ib$ (62)

so that

${\left(1+\mu -a\right)}^{2}+{b}^{2}=-\left(\alpha \delta +\beta \sigma \right)$ (63)

and also that $\mathrm{det}U=-b\left(\alpha \delta +\beta \sigma \right)$. So if

$\frac{\mathrm{sin}\left({b}_{1}\left(t-T\right)\right)}{-b}>0$ (64)

then

${u}^{*}\left(t\right)={g}_{I}^{0}$ (65)

is optimal and if the reverse inequality holds then

${u}^{*}\left(t\right)=0$ (66)

is optimal. We shall now consider the case where all eigenvalues ${\lambda }_{+},{\lambda }_{-},\stackrel{˜}{\lambda }$ are real, positive and distinct and that $\mu ={\mu }_{F}={\mu }_{I}$. Here we let

$Y\left(x,u\right)=\stackrel{˜}{\Lambda }x+d+{D}^{-1}{e}_{3}u$ (67)

where

$\stackrel{˜}{\Lambda }=\text{diag}\left(\mathrm{ln}{\lambda }_{+},\mathrm{ln}{\lambda }_{-},\mathrm{ln}\left(1+\mu \right)\right)$ (68)

Also define d in

$\frac{{\lambda }_{+}-1}{\mathrm{ln}{\lambda }_{+}}{d}_{1}={\left({D}^{-1}g\right)}_{1}$ (69)

$\frac{{\lambda }_{-}-1}{\mathrm{ln}{\lambda }_{+}}{d}_{2}={\left({D}^{-1}g\right)}_{2}$ (70)

$\frac{\mu }{\mathrm{ln}\left(1+\mu \right)}{d}_{3}={\left({D}^{-1}g\right)}_{3}$ (71)

$d={\left({d}_{1},{d}_{2},{d}_{3}\right)}^{\text{T}}\in {ℝ}^{3}$. Now define the Hamiltonian

$H\left(x,p,u,t\right)={p}^{\text{T}}Y\left(x,u\right)$ (72)

Then we get the adjoint equations

${{p}^{\prime }}_{1}=-\frac{\partial H}{\partial {x}_{1}}=-\mathrm{ln}{\lambda }_{+}{p}_{1}$ (73)

${{p}^{\prime }}_{2}=-\frac{\partial H}{\partial {x}_{2}}=-\mathrm{ln}{\lambda }_{-}{p}_{2}$ (74)

${{p}^{\prime }}_{3}=-\frac{\partial H}{\partial {x}_{3}}=-\mathrm{ln}\left(1+\mu \right){p}_{3}$ (75)

see . Now define

${S}_{1}\left(x\right)=\left(1+\mu -{\lambda }_{+}\right){x}_{1}+\left(1+\mu -{\lambda }_{-}\right){x}_{2}$ (76)

Then we have the transversality conditions

$p\left(T\right)=\frac{\partial {S}_{1}}{\partial x}=\left(\begin{array}{c}1+\mu -{\lambda }_{+}\\ 1+\mu -{\lambda }_{-}\\ 0\end{array}\right)$ (77)

by  and see below. Now observe, that

$\mathrm{det}\left(D\right)=\left({\lambda }_{+}-{\lambda }_{-}\right)\left(\alpha \delta +\beta \sigma \right)$ (78)

Because

$1+\mu -{\lambda }_{±}=\frac{\mu -\gamma }{2}\mp \frac{1}{2}\sqrt{{\left(\mu -\gamma \right)}^{2}+4\left(\alpha \delta +\beta \sigma \right)}$ (79)

we find

$\left(1+\mu -{\lambda }_{+}\right)\left(1+\mu -{\lambda }_{-}\right)=-\left(\alpha \delta +\beta \sigma \right)$ (80)

We shall need

${D}^{-1}=\left(\begin{array}{ccc}-\alpha \delta -\beta \sigma & -\alpha \left(1+\mu -{\lambda }_{-}\right)& -\beta \left(1+\mu -{\lambda }_{-}\right)\\ \alpha \delta +\beta \sigma & \alpha \left(1+\mu -{\lambda }_{+}\right)& \beta \left(1+\mu -{\lambda }_{+}\right)\\ 0& -\sigma \left({\lambda }_{+}-{\lambda }_{-}\right)& \delta \left({\lambda }_{+}-{\lambda }_{-}\right)\end{array}\right)\frac{1}{\mathrm{det}D}$ (81)

We have

$\left(\begin{array}{c}C\left(T\right)\\ 0\\ 0\end{array}\right)=D\left(\begin{array}{c}{x}_{1}\left(T\right)\\ {x}_{2}\left(T\right)\\ {x}_{3}\left(T\right)\end{array}\right)$ (82)

hence the definition of ${S}_{1}$. Now

$p\left(T\right)=\frac{\partial {S}_{1}}{\partial x}$ (83)

thus

${p}_{1}\left(T\right)={{p}_{1}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}{\lambda }_{+}t\right)|}_{t=T}={p}_{1}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}{\lambda }_{+}T\right)=1+\mu -{\lambda }_{+}$ (84)

${p}_{2}\left(T\right)={{p}_{2}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}{\lambda }_{-}t\right)|}_{t=T}={p}_{2}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}{\lambda }_{-}T\right)=1+\mu -{\lambda }_{-}$ (85)

${p}_{3}\left(T\right)={{p}_{3}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}\left(1+\mu \right)t\right)|}_{t=T}={p}_{3}\left(0\right)\mathrm{exp}\left(-\mathrm{ln}\left(1+\mu \right)T\right)=0$ (86)

So

${p}_{1}\left(t\right)=\left(1+\mu -{\lambda }_{+}\right)\mathrm{exp}\left(\mathrm{ln}{\lambda }_{+}\left(T-t\right)\right)$ (87)

${p}_{2}\left(t\right)=\left(1+\mu -{\lambda }_{-}\right)\mathrm{exp}\left(\mathrm{ln}{\lambda }_{-}\left(T-t\right)\right)$ (88)

${p}_{3}\left(t\right)=0$ (89)

If

$H\left({x}^{*}\left(t\right),p\left(t\right),{u}^{*}\left(t\right),t\right)\le H\left({x}^{*}\left(t\right),p\left(t\right),u,t\right)$ (90)

for all $u\in \mathbb{U}$,$\left({x}^{*}\left(t\right),{u}^{*}\left(t\right)\right)$ is optimal. And this is equivalent to

$p{\left(t\right)}^{\text{T}}{D}^{-1}{e}_{3}{u}^{*}\left(t\right)\le p{\left(t\right)}^{\text{T}}{D}^{-1}{e}_{3}u$ (91)

which again is equivalent to

$\left(1+\mu -{\lambda }_{+}\right)\mathrm{exp}\left(\mathrm{ln}{\lambda }_{+}\left(T-t\right)\right)\left(-\beta \right)\left(1+\mu -{\lambda }_{-}\right)\frac{{u}^{*}\left(t\right)}{\mathrm{det}\left(D\right)}$ (92)

$+\left(1+\mu -{\lambda }_{-}\right)\mathrm{exp}\left(\mathrm{ln}{\lambda }_{-}\left(T-t\right)\right)\beta \left(1+\mu -{\lambda }_{+}\right)\frac{{u}^{*}\left(t\right)}{\mathrm{det}\left(D\right)}$ (93)

$={u}^{*}\left(t\right)\beta \frac{-\alpha \delta -\beta \sigma }{\left({\lambda }_{+}-{\lambda }_{-}\right)\left(\alpha \delta +\beta \sigma \right)}\left(-\mathrm{exp}\left(\mathrm{ln}{\lambda }_{+}\left(T-t\right)\right)+\mathrm{exp}\left(\mathrm{ln}{\lambda }_{-}\left(T-t\right)\right)\right)$ (94)

$\le u\beta \left(\mathrm{exp}\left(\mathrm{ln}{\lambda }_{+}\left(T-t\right)\right)-\mathrm{exp}\left(\mathrm{ln}{\lambda }_{-}\left(T-t\right)\right)\right)\frac{1}{{\lambda }_{+}-{\lambda }_{-}}$ (95)

for all $u\in \mathbb{U}$. It follows that

${u}^{*}\left(t\right)={g}_{I}^{0}$ (96)

is optimal. We have the two Hamiltonians

${H}^{Y}\left(x,p,u,t\right)={p}^{\text{T}}Y\left(x,u\right)=〈p,Y\left(x,u\right)〉$ (97)

${H}^{X}\left(y,q,v,t\right)={q}^{\text{T}}X\left(y,v\right)=〈q,X\left(y,v\right)〉$ (98)

where $〈\cdot ,\cdot 〉$ denotes the canonical inner product. It follows that when the eigenvalues are $a±ib,1+\mu >0,a>0,b\ne 0$

$y\left(t\right)=Ux\left(t\right)$ (99)

$p\left(t\right)={U}^{\text{T}}q\left(t\right)$ (100)

if

$y\left(0\right)=Ux\left(0\right)$ (101)

$p\left(T\right)={U}^{\text{T}}q\left(T\right)$ (102)

where $y\left(t\right)$ is an integral curve of X and $x\left(t\right)$ an integral curve of Y. Because

${q}^{\prime }=-\frac{\partial {H}^{X}}{\partial y}=-{\left({U}^{-1}\right)}^{\text{T}}{\stackrel{˜}{B}}^{\text{T}}{U}^{\text{T}}q$ (103)

hence

${U}^{\text{T}}{q}^{\prime }=-{\stackrel{˜}{B}}^{\text{T}}{U}^{\text{T}}q$ (104)

and

${p}^{\prime }=-{\stackrel{˜}{B}}^{\text{T}}p$ (105)

We now get, that

${H}^{X}\left(y\left(t\right),q\left(t\right),v\left(t\right),t\right)=〈q\left(t\right),X\left(y\left(t\right),v\left(t\right)\right)〉$ (106)

$=〈q\left(t\right),UY\left({U}^{-1}y\left(t\right),v\left(t\right)\right)〉$ (107)

$=〈{U}^{\text{T}}\left(q\left(t\right)\right),Y\left({U}^{-1}\left(y\left(t\right)\right),v\left(t\right)\right)〉$ (108)

$=〈p\left(t\right),Y\left(x\left(t\right),v\left(t\right)\right)〉$ (109)

$={H}^{Y}\left(x\left(t\right),p\left(t\right),v\left(t\right),t\right)$ (110)

So

${H}^{X}\left(y\left(t\right),q\left(t\right),v\left(t\right),t\right)={H}^{Y}\left(x\left(t\right),p\left(t\right),v\left(t\right),t\right)$ (111)

When the eigenvalues of A are real, distinct and positive

$y\left(t\right)=Dx\left(t\right)$ (112)

$p\left(t\right)={D}^{\text{T}}q\left(t\right)$ (113)

if

$y\left(0\right)=Dx\left(0\right)$ (114)

$p\left(T\right)={D}^{\text{T}}q\left(T\right)$ (115)

We also have

${q}^{\prime }=-{\left({D}^{-1}\right)}^{\text{T}}{\Lambda }^{\text{T}}{D}^{\text{T}}q$ (116)

and

${\left({D}^{\text{T}}q\right)}^{\prime }=-{\Lambda }^{\text{T}}{D}^{\text{T}}q$ (117)

thus

${p}^{\prime }=-{\Lambda }^{\text{T}}p$ (118)

Hence

${H}^{X}\left(y\left(t\right),q\left(t\right),v\left(t\right),t\right)={H}^{Y}\left(x\left(t\right),p\left(t\right),v\left(t\right),t\right)$ (119)

We need the following theorem, which is well known.

Theorem 1 Now The following statements about a C2 function from

$f:{ℝ}^{n}\to ℝ$ (120)

to the reals, where $n$ is a positive integer, $n\in ℕ$ are equivalent:

(i) $f\left(\lambda x+\left(1-\lambda \right)y\right)\le \lambda f\left(x\right)+\left(1-\lambda \right)f\left(y\right)$ where $\lambda \in \left[0,1\right]$ ;

(ii) $f\left(y\right)\ge f\left(x\right)+\frac{\partial f}{\partial x}\left(x\right)\left(y-x\right)$ ;

(iii) ${\sum }_{i,j=1,\cdots ,n}\text{ }{\gamma }_{i}\frac{{\partial }^{2}f}{\partial {x}_{i}\partial {x}_{j}}{\gamma }_{j}\ge 0$

where $x,y,\gamma \in {ℝ}^{n}$.

$\left(x\left(t\right),u\left(t\right)\right)$ is admissible by definition if $0\le u\left(t\right)\le {g}_{I}^{0}$ and

${x}^{\prime }\left(t\right)=Y\left(x\left(t\right),u\left(t\right)\right),\text{\hspace{0.17em}}\text{\hspace{0.17em}}x\left(0\right)\text{=}{x}_{0}\in {ℝ}^{3}$ (121)

To see that

${S}_{1}\left({x}^{*}\left(T\right)\right)\le {S}_{1}\left(x\left(T\right)\right)$ (122)

for $\left({x}^{*}\left(t\right),{u}^{*}\left(t\right)\right)$ optimal candidate and $\left(x\left(t\right),u\left(t\right)\right)$ admissible argue as in 

$\Delta ={S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)$ (123)

$={S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)$ (124)

$\begin{array}{l}+{\int }_{0}^{T}\left({H}^{Y}\left({x}^{*}\left(t\right),p\left(t\right),{u}^{*}\left(t\right),t\right)-p{\left(t\right)}^{\text{T}}{x}^{*}{}^{\prime }\left(t\right)\right)\\ -\left({H}^{Y}\left(x\left(t\right),p\left(t\right),u\left(t\right),t\right)-p{\left(t\right)}^{\text{T}}{x}^{\prime }\left(t\right)\right)\text{d}t\end{array}$ (125)

We have the following inequality, which follows from theorem 1.

$\begin{array}{l}{H}^{Y}\left({x}^{*},p,{u}^{*},t\right)-{H}^{Y}\left(x,p,u,t\right)\\ \le \frac{\partial {H}^{Y}}{\partial x}\left({x}^{*},p,{u}^{*},t\right)\left({x}^{*}-x\right)+\frac{\partial {H}^{Y}}{\partial u}\left({x}^{*},p,{u}^{*},t\right)\left({u}^{*}-u\right)\end{array}$ (126)

We also have

${p}^{\prime }=-\frac{\partial {H}^{Y}}{\partial x}\left({x}^{*},p,{u}^{*},t\right)$ (127)

We can now estimate

$\Delta \le {S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)+{\int }_{0}^{T}{{p}^{\prime }}^{\text{T}}\left(x-{x}^{*}\right)+{p}^{\text{T}}\left({x}^{\prime }-{x}^{*}{}^{\prime }\right)\text{d}t$ (128)

$\begin{array}{l}\text{ }+{\int }_{0}^{T}\frac{\partial {H}^{Y}}{\partial u}\left({x}^{*},p,{u}^{*},t\right)\left({u}^{*}-u\right)\text{d}t\\ \le {\int }_{0}^{T}\frac{\text{d}}{\text{d}t}\left({p}^{\text{T}}\left(x-{x}^{*}\right)\right)\left(t\right)\text{d}t+{S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)\end{array}$ (129)

$=p{\left(T\right)}^{\text{T}}\left(x\left(T\right)-{x}^{*}\left(T\right)\right)+{S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)$ (130)

$=\frac{\partial {S}_{1}}{\partial x}\left({x}^{*}\left(T\right)\right)\left(x\left(T\right)-{x}^{*}\left(T\right)\right)+{S}_{1}\left({x}^{*}\left(T\right)\right)-{S}_{1}\left(x\left(T\right)\right)$ (131)

$\le 0$ (132)

because ${S}_{1}$ is convex, by theorem 1. We have used, that we have arranged, that

$\frac{\partial H}{\partial u}\left({x}^{*},p,{u}^{*},t\right)\left({u}^{*}-u\right)\le 0$ (133)

for all $u\in \mathbb{U}$, by the mean value theorem.

We have optimality.

3. Optimal Control of T

In this section we consider the problem: minimize $C\left(N\right)$ subject to

${y}_{k+1}=A{y}_{k}+B{u}_{k}+g=f\left({y}_{k},{u}_{k}\right)$ (134)

$k=1,\cdots ,N-1,N\in ℕ,N\ge 2,u\left(k\right)\in \mathbb{U}=\left[0,{g}_{I}^{0}\right],{g}_{i}^{0}>0,\mu ={\mu }_{F}={\mu }_{I}$ where A is as in the introduction.

Also

${y}_{k}=\left(C\left(k\right),GF\left(k\right),GI\left(k\right)\right),g\in {ℝ}^{3}$ (135)

Here

$B={\left(0,0,1\right)}^{\text{T}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}g={\left({g}_{C},{g}_{F},{g}_{I}\right)}^{\text{T}}$ (136)

Assume (i) the eigenvalues of A are real and distinct.

In the Discrete Pontryagin Minimal Principle applied to T you define the Hamiltonian by (138) and then you find

$\frac{\partial H}{\partial {u}_{k}}\left({x}_{k}^{*},{u}_{k}^{*},{\lambda }_{k}\right)\left({u}_{k}^{*}-{u}_{k}\right)$ (137)

and minimize it to find the optimal control ${u}_{k}^{*}$. It is optimal due to computations (157) to (163) below.

Define then the Hamiltonian

$H\left({x}_{k},{u}_{k},{\lambda }_{k}\right)={\lambda }_{k}^{\text{T}}\left(\stackrel{˜}{A}{x}_{k}+\stackrel{^}{B}{u}_{k}+{D}^{-1}g\right)$ (138)

where ${\lambda }_{k}\in {ℝ}^{3}$ and (139)

$\stackrel{^}{B}={D}^{-1}B$ (140)

$\stackrel{˜}{A}={D}^{-1}AD$ (141)

Then we have the adjoint equation

${\lambda }_{k-1}=\frac{\partial H}{\partial {x}_{k}}\left({x}_{k},{u}_{k},{\lambda }_{k}\right)={\stackrel{˜}{A}}^{\text{T}}{\lambda }_{k}$ (142)

Inductively

${\lambda }_{N-k-1}={\left({\stackrel{˜}{A}}^{\text{T}}\right)}^{k}{\lambda }_{N-1}$ (143)

In particular

${\lambda }_{0}={\left({\stackrel{˜}{A}}^{\text{T}}\right)}^{N-1}{\lambda }_{N-1}$ (144)

For $k=0$ we have

$\begin{array}{l}H\left({x}_{0}^{*},{u}_{0}^{*},{\lambda }_{0}\right)={\left({\lambda }_{0}\right)}^{\text{T}}\left(\stackrel{˜}{A}{x}_{0}^{*}+\stackrel{^}{B}{u}_{0}^{*}+{D}^{-1}g\right)\\ \le H\left({x}_{0}^{*},u,{\lambda }_{0}\right)={\left({\lambda }_{0}\right)}^{\text{T}}\left(\stackrel{˜}{A}{x}_{0}^{*}+\stackrel{^}{B}u+{D}^{-1}g\right)\end{array}$ (145)

which is equivalent to

${\left({\lambda }_{0}\right)}^{\text{T}}\stackrel{^}{B}{u}_{0}^{*}\le {\left({\lambda }_{0}\right)}^{\text{T}}\stackrel{^}{B}u$ (146)

Define

${S}_{1}\left(x\right)=F\left(x\right)=\left(1+\mu -{\lambda }_{+}\right){x}_{1}+\left(1+\mu -{\lambda }_{-}\right){x}_{2}$ (147)

Now

${\lambda }_{0}={\left({\stackrel{˜}{A}}^{\text{T}}\right)}^{N-1}{\lambda }_{N-1}$ (148)

where

${\lambda }_{N-1}=\frac{\partial F}{\partial x}\left({x}_{N}\right)$ (149)

Thus

${\left({\lambda }_{0}\right)}^{\text{T}}\stackrel{^}{B}{u}_{0}^{*}={\left({\left({\stackrel{˜}{A}}^{\text{T}}\right)}^{N-1}{\lambda }_{N-1}\right)}^{\text{T}}\stackrel{^}{B}{u}_{0}^{*}={\left({\lambda }_{N-1}\right)}^{\text{T}}{\stackrel{˜}{A}}^{N-1}\stackrel{^}{B}{u}_{0}^{*}={\frac{\partial F}{\partial x}}^{\text{T}}\left({x}_{N}\right){\stackrel{˜}{A}}^{N-1}\stackrel{^}{B}{u}_{0}^{*}$ (150)

Assume that (i) holds and ${\lambda }_{+},{\lambda }_{-},1+\mu \in ℝ$ are distinct, when $\mu ={\mu }_{F}={\mu }_{I}$. Then

$\stackrel{^}{B}=\left(\begin{array}{c}{D}_{31}\\ {D}_{32}\\ {D}_{33}\end{array}\right)\frac{1}{\mathrm{det}\left(D\right)}=\left(\begin{array}{c}-\beta \left(1+\mu -{\lambda }_{-}\right)\\ \beta \left(1+\mu -{\lambda }_{+}\right)\\ \delta \left({\lambda }_{+}-{\lambda }_{-}\right)\end{array}\right)\frac{1}{\mathrm{det}\left(D\right)}={D}^{-1}{e}_{3}$ (151)

So now we get, that (146) amounts to

$\left(1+\mu -{\lambda }_{+},1-\mu -{\lambda }_{-},0\right)\left(\begin{array}{ccc}{\lambda }_{+}^{N-1}& 0& 0\\ 0& {\lambda }_{-}^{N-1}& 0\\ 0& 0& {\left(1+\mu \right)}^{N-1}\end{array}\right)\stackrel{^}{B}$ (152)

$=\frac{-\beta \left(-\alpha \delta -\beta \sigma \right){\lambda }_{+}^{N-1}}{\mathrm{det}\left(D\right)}+\frac{\beta \left(-\alpha \delta -\beta \sigma \right){\lambda }_{-}^{N-1}}{\mathrm{det}\left(D\right)}$ (153)

$=\frac{\alpha \delta +\beta \sigma }{\left({\lambda }_{+}-{\lambda }_{-}\right)\left(\alpha \delta +\beta \sigma \right)}\beta \left({\lambda }_{+}^{N-1}-{\lambda }_{-}^{N-1}\right)$ (154)

$=\frac{\beta }{{\lambda }_{+}-{\lambda }_{-}}\left({\lambda }_{+}^{N-1}-{\lambda }_{-}^{N-1}\right)\le 0$ (155)

Similarly for $k=1,\cdots ,N-2$

${\lambda }_{k}^{\text{T}}\stackrel{^}{B}{u}_{k}^{*}=\left(1+\mu -{\lambda }_{+},1+\mu -{\lambda }_{-},0\right)\left(\begin{array}{ccc}{\lambda }_{+}^{N-k-1}& 0& 0\\ 0& {\lambda }_{-}^{N-k-1}& 0\\ 0& 0& {\left(1+\mu \right)}^{N-k-1}\end{array}\right)\stackrel{^}{B}{u}_{k}^{*}<0$ (156)

For $k=N-1$ we have, that $\left(1,0,0\right)\cdot \stackrel{^}{B}=0$. This means that maximal chemo therapy is optimal. Because similar to (128) to (132), we get

$\Delta ={S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)$ (157)

$={S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)+\underset{k=0}{\overset{N-1}{\sum }}\left(H\left({x}_{k}^{*},{u}_{k}^{*},{\lambda }_{k}\right)-H\left({x}_{k},{u}_{k},{\lambda }_{k}\right)-{\lambda }_{k}^{\text{T}}{x}_{k+1}^{*}+{\lambda }_{k}^{\text{T}}{x}_{k+1}\right)$ (158)

$\le {S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)+\underset{k=0}{\overset{N-1}{\sum }}\frac{\partial H}{\partial {x}_{k}}\left({x}_{k}^{*},{u}_{k}^{*},{\lambda }_{k}\right)\left({x}_{k}^{*}-{x}_{k}\right)$ (159)

$\text{+}\frac{\partial H}{\partial {u}_{k}}\left({x}_{k}^{*},{u}_{k}^{*},{\lambda }_{k}\right)\left({u}_{k}^{*}-{u}_{k}\right)-{\lambda }_{k}^{\text{T}}{x}_{k+1}^{*}+{\lambda }_{k}^{\text{T}}{x}_{k+1}$ (160)

$\le {S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)+\underset{k=0}{\overset{N-1}{\sum }}\text{ }{\lambda }_{k-1}^{\text{T}}\left({x}_{k}^{*}-{x}_{k}\right)-{\lambda }_{k}^{\text{T}}{x}_{k+1}^{*}+{\lambda }_{k}^{\text{T}}{x}_{k+1}$ (161)

$={S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)-{\lambda }_{N-1}^{\text{T}}\left({x}_{N}^{*}-{x}_{N}\right)$ (162)

$={S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)-\frac{\partial {S}_{1}}{\partial x}\left({x}_{N}^{*}\right)\left({x}_{N}^{*}-{x}_{N}\right)\le 0$ (163)

by theorem 1 and since ${S}_{1}$ is convex.

Then we have

${C}_{*}\left(N\right)-C\left(N\right)={S}_{1}\left({x}_{N}^{*}\right)-{S}_{1}\left({x}_{N}\right)\le 0$ (164)

As above we get

$\stackrel{˜}{H}\left({y}_{k},{u}_{k},{\zeta }_{k}\right)={\zeta }_{k}^{\text{T}}\left(A{y}_{k}+g+{e}_{3}{v}_{k}\right)$ (165)

and

${y}_{k+1}=A{y}_{k}+g+{e}_{3}{v}_{k}$ (166)

Also

${\zeta }_{k-1}={A}^{\text{T}}{\zeta }_{k}$ (167)

and

${\lambda }_{k-1}={\stackrel{˜}{A}}^{\text{T}}{\lambda }_{k}$ (168)

Thus

${\lambda }_{k-1}={D}^{\text{T}}{\zeta }_{k-1}$ (169)

When

${y}_{0}=D{x}_{0}$ (170)

then

${y}_{k}=D{x}_{k}$ (171)

Hence

$\stackrel{˜}{H}\left({y}_{k},{v}_{k},{\zeta }_{k}\right)=H\left({x}_{k},{v}_{k},{\lambda }_{k}\right)$ (172)

Now consider the case where there are imaginary eigenvalues $a±ib,a>0,b\ne 0,1+\mu >0$ for A. We need the following well known formulas for

${\left(\begin{array}{cc}a& b\\ -b& a\end{array}\right)}^{2p+1}=\left(\begin{array}{cc}{A}_{p}& {B}_{p}\\ -{B}_{p}& {A}_{p}\end{array}\right)$ (173)

which are well known, where $p\in {ℕ}_{0}$ and

${A}_{p}=\underset{q=0}{\overset{p}{\sum }}\left(\begin{array}{c}2p+1\\ 2q\end{array}\right){\left(-1\right)}^{q}{b}^{2q}{a}^{2p+1-2q}$ (174)

${B}_{p}=\underset{q=0}{\overset{p}{\sum }}\left(\begin{array}{c}2p+1\\ 2q+1\end{array}\right){\left(-1\right)}^{q}{b}^{2q+1}{a}^{2p-2q}$ (175)

and $p\in ℕ$

${\left(\begin{array}{cc}a& b\\ -b& a\end{array}\right)}^{2p}=\left(\begin{array}{cc}{C}_{p}& {D}_{p}\\ -{D}_{p}& {C}_{p}\end{array}\right)$ (176)

where

${C}_{p}=\underset{q=0}{\overset{p}{\sum }}\left(\begin{array}{c}2p\\ 2q\end{array}\right){\left(-1\right)}^{q}{b}^{2q}{a}^{2p-2q}$ (177)

${D}_{p}=\underset{q=0}{\overset{p-1}{\sum }}\left(\begin{array}{c}2p\\ 2q+1\end{array}\right){\left(-1\right)}^{q}{b}^{2q+1}{a}^{2p-1-2q}$ (178)

You can prove them by induction. We have

$\begin{array}{l}H\left({x}_{k},{v}_{k},{\lambda }_{k}\right)=〈{\lambda }_{k},\stackrel{˜}{B}{x}_{k}+{U}^{-1}{e}_{3}{v}_{k}+{U}^{-1}g〉\\ =\stackrel{˜}{H}\left({y}_{k},{v}_{k},{\zeta }_{k}\right)=〈{\zeta }_{k},A{y}_{k}+g+{e}_{3}{v}_{k}〉\end{array}$ (179)

where

$\stackrel{˜}{B}=\left(\begin{array}{ccc}a& b& 0\\ -b& a& 0\\ 0& 0& 1+\mu \end{array}\right)$ (180)

As above

${\lambda }_{k}={\left({\stackrel{˜}{B}}^{\text{T}}\right)}^{N-k-1}{\lambda }_{N-1}={\left({\stackrel{˜}{B}}^{\text{T}}\right)}^{N-k-1}\frac{\partial F}{\partial x}$ (181)

where

${S}_{1}\left(x\right)=F\left(x\right)=\left(1+\mu -a\right){x}_{1}-b{x}_{2}$ (182)

So

${\lambda }_{k}^{\text{T}}\stackrel{^}{B}=\left(1+\mu -a,-b,0\right){\stackrel{˜}{B}}^{N-k-1}{\left(\beta b，\beta \left(1+\mu -a\right)，*\right)}^{\text{T}}\frac{1}{\mathrm{det}\left(U\right)}$ (183)

For $k\ge 0$ we find when $N-k-1=2p+1,k\ne N-1$

$\frac{\partial H}{\partial {v}_{k}}=\left(1+\mu -a,-b\right)\left(\begin{array}{cc}{A}_{p}& {B}_{p}\\ -{B}_{p}& {A}_{p}\end{array}\right)\beta \left(\begin{array}{c}b\\ 1+\mu -a\end{array}\right)\frac{1}{\mathrm{det}\left(U\right)}$ (184)

$=\beta \left(\left(1+\mu -a\right)\left({A}_{p}b+\left(1+\mu -a\right){B}_{p}\right)-b\left(-{B}_{p}b+\left(1+\mu -a\right){A}_{p}\right)\right)\frac{1}{\text{det}\left(U\right)}$ (185)

$=\beta \left({\left(1+\mu -a\right)}^{2}+{b}^{2}\right){B}_{p}\frac{1}{\mathrm{det}\left(U\right)}$ (186)

$=\beta \left(-\alpha \delta -\beta \sigma \right){B}_{p}\frac{1}{\mathrm{det}\left(U\right)}$ (187)

$=\frac{\beta }{b}{B}_{p}$ (188)

Here $\mathrm{det}\left(U\right)=-b\left(\alpha \delta +\beta \sigma \right)$. When $N-k-1=2p,k\ne N-1$

$\frac{\partial H}{\partial {v}_{k}}=\left(1+\mu -a,-b\right)\left(\begin{array}{cc}{C}_{p}& {D}_{p}\\ -{D}_{p}& {C}_{p}\end{array}\right)\beta \left(\begin{array}{c}b\\ 1+\mu -a\end{array}\right)\frac{1}{\text{det}\left(U\right)}=\frac{\beta }{b}{D}_{p}$ (189)

If

$\frac{\beta }{b}{B}_{p}<0$ (190)

let

${v}_{k}^{*}={g}_{I}^{0}$ (191)

and

${v}_{k}^{*}=0$ (192)

if the reverse inequality holds. If

$\frac{\beta }{b}{D}_{p}<0$ (193)

let

${v}_{k}^{*}={g}_{I}^{0}$ (194)

and

${v}_{k}^{*}=0$ (195)

if the reverse inequality holds. Then $\left({x}_{k}^{*},{v}_{k}^{*}\right)$ is optimal, by (157) to (163). For $k=N-1$ we have $\left(1,0,0\right)\cdot B=0$.

4. A Counter Example

We shall now present a counter example to optimality of maximal chemo therapy when the eigenvalues are real and ${\mu }_{F}\ne {\mu }_{I}$ for the model

$T\left(y\right)=Ay+g$ (196)

of the introduction.

Remember the definition of the discriminant $\Delta$ of a cubic polynomial

$p\left(\lambda \right)={\lambda }^{3}+{a}_{1}{\lambda }^{2}+{a}_{2}\lambda +{a}_{3}$ (197)

Namely

$-108\Delta ={a}_{1}^{2}{a}_{2}^{2}-27{a}_{3}^{2}-4{a}_{2}^{3}-4{a}_{1}^{3}{a}_{3}+18{a}_{1}{a}_{2}{a}_{3}$ (198)

Notice that the degree of $1+\gamma$ is four in ${a}_{1}^{2}{a}_{2}^{2}$ and $-4{a}_{1}^{3}{a}_{3}$ while it is only two or three in $-27{a}_{3}^{2},-4{a}_{2}^{3},18{a}_{1}{a}_{2}{a}_{3}$. We thus get

$-108\Delta ={a}_{1}^{2}{a}_{2}^{2}-4{a}_{1}^{3}{a}_{3}+\text{lowerordertermsin}\left(1+\gamma \right)$ (199)

$={\left(3+\gamma +{\mu }_{F}+{\mu }_{I}\right)}^{2}\left(\alpha \delta +\beta \sigma \right)-\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)-\left(1+\gamma \right){\left(2+{\mu }_{F}+{\mu }_{I}\right)}^{2}$ (200)

$\text{ }+4{\left(3+\gamma +{\mu }_{F}+{\mu }_{I}\right)}^{3}\left(\alpha \delta \left(1+{\mu }_{I}\right)+\beta \sigma \left(1+{\mu }_{F}\right)-\left(1+\gamma \right)\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)\right)$ (201)

$\text{ }+\text{lowerordertermsin}\left(1+\gamma \right)$ (202)

which becomes

$\left({\left(1+\gamma \right)}^{2}+{\left(1+{\mu }_{F}\right)}^{2}+{\left(1+{\mu }_{I}\right)}^{2}+2\left(1+\gamma \right)\left(1+{\mu }_{F}\right)+2\left(1+\gamma \right)\left(1+{\mu }_{I}\right)$ (203)

$\text{ }+2\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)\right){\left(\alpha \delta +\beta \sigma -\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)-\left(1+\gamma \right)\left(2+{\mu }_{F}+{\mu }_{I}\right)\right)}^{2}$ (204)

$\begin{array}{l}\text{ }+4\left({\left(1+\gamma \right)}^{3}+\text{lowerordertermsin}\left(1+\gamma \right)\right)\left(\alpha \delta \left(1+{\mu }_{I}\right)+\beta \sigma \left(1+{\mu }_{F}\right)\\ \text{ }-\left(1+\gamma \right)\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)\right)\end{array}$ 205)

$\begin{array}{l}={\left(1+\gamma \right)}^{2}{\left(1+\gamma \right)}^{2}{\left(1+{\mu }_{F}+1+{\mu }_{I}\right)}^{2}-4{\left(1+\gamma \right)}^{4}\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{ }\text{ }+\text{lowerordertermsin}\left(1+\gamma \right)\end{array}$ (206)

and this is

$\begin{array}{l}{\left(1+\gamma \right)}^{4}\left({\left(1+{\mu }_{F}\right)}^{2}+{\left(1+{\mu }_{I}\right)}^{2}+2\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)-4\left(1+{\mu }_{F}\right)\left(1+{\mu }_{I}\right)\right)\\ \text{ }+\text{lowerordertermsin}\left(1+\gamma \right)\end{array}$ (207)

$={\left(1+\gamma \right)}^{4}{\left({\mu }_{F}-{\mu }_{I}\right)}^{2}+\text{lowerordertermsin}\left(1+\gamma \right)$ (208)

It follows that for $-\gamma$ large

$\Delta <0$ (209)

so that there are three real distinct roots, see Uspensky . But

${A}_{13}^{2}=\beta \left(1+\gamma \right)+\beta \left(1+{\mu }_{I}\right)>0$ (210)

when $-\gamma$ is large. So maximal chemo therapy is not optimal for $N=3$. In fact

${x}_{3}={A}^{3}{x}_{0}+{A}^{2}\left(B{u}_{0}+g\right)+A\left(B{u}_{1}+g\right)+B{u}_{2}+g$ (211)

So ${u}_{0}=0,{u}_{1}={g}_{I}^{0}$ gives an optimal trajectory.

5. Optimality of T When ${\mu }_{F}\ne {\mu }_{I}$

Consider the model T from the introduction

$T\left(y\right)=Ay+g+u{e}_{3}$ (212)

where

$A=\left(\begin{array}{ccc}1+\gamma & \alpha & \beta \\ \delta & 1+{\mu }_{F}& 0\\ \sigma & 0& 1+{\mu }_{I}\end{array}\right)$ (213)

and

$y={\left(C,GF,GI\right)}^{\text{T}},\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}g={\left({g}_{C},{g}_{F},{g}_{I}\right)}^{\text{T}}\in {ℝ}^{3}$ (214)

Assume (i) the eigenvalues ${\lambda }_{+},{\lambda }_{-},\stackrel{˜}{\lambda }=1+\mu$ of A are real and distinct, when $\mu ={\mu }_{F}={\mu }_{I}$.

Theorem 2 There exists a Euclidean open ball ${B}_{\rho }\left(\mu ,\mu \right),\rho >0$ in ${ℝ}^{2}$, such that for $\left({\mu }_{F},{\mu }_{I}\right)\in {B}_{\rho }\left(\mu ,\mu \right)$ maximal chemo therapy is optimal.

Also let

$D=\left(\begin{array}{ccc}\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)& \left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)& -\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\\ -\delta \left(1+{\mu }_{I}-{\lambda }_{+}\right)& -\delta \left(1+{\mu }_{I}-{\lambda }_{-}\right)& \delta \beta \\ -\sigma \left(1+{\mu }_{F}-{\lambda }_{+}\right)& -\sigma \left(1+{\mu }_{F}-{\lambda }_{-}\right)& -\alpha \delta +\left(1+\gamma -\stackrel{˜}{\lambda }\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\end{array}\right)$ (215)

be a matrix with column eigenvectors to the eigenvalues ${\lambda }_{+},{\lambda }_{-},\stackrel{˜}{\lambda }$ of A. We have ( $\delta \ne 0$ )

${\text{det}\left(D\right)|}_{{\mu }_{F}={\mu }_{I}=\mu }=\left(1+\mu -{\lambda }_{+}\right)\left(1+\mu -{\lambda }_{-}\right)|\begin{array}{ccc}1+\mu -{\lambda }_{+}& 1+\mu -{\lambda }_{-}& 0\\ -\delta & -\delta & \beta \delta \\ -\sigma & -\sigma & -\alpha \delta \end{array}|$ (216)

$=\left(1+\mu -{\lambda }_{+}\right)\left({\delta }^{2}\alpha +\beta \sigma \delta \right)\left(-\alpha \delta -\sigma \beta \right)$ (217)

$-\left(1+\mu -{\lambda }_{-}\right)\left({\delta }^{2}\alpha +\beta \sigma \delta \right)\left(-\alpha \delta -\sigma \beta \right)$ (218)

$=\left({\lambda }_{+}-{\lambda }_{-}\right){\left(\alpha \delta +\beta \sigma \right)}^{2}\delta$ (219)

Define the Hamiltonian

$H\left({x}_{k},{u}_{k},{\lambda }_{k}\right)={\lambda }_{k}^{\text{T}}{f}_{k}\left({x}_{k},{u}_{k}\right)={\lambda }_{k}^{\text{T}}\left(A{x}_{k}+g+B{u}_{k}\right)$ (220)

$k=0,\cdots ,N-1$ and

$\begin{array}{c}F\left(x\right)=\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right){x}_{1}+\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right){x}_{2}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}-\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right){x}_{3}\end{array}$ (221)

We have

${\lambda }_{N-k-1}={\left({A}^{\text{T}}\right)}^{k}{\lambda }_{N-1}$ (222)

In particular

${\lambda }_{0}={\left({A}^{\text{T}}\right)}^{N-1}{\lambda }_{N-1}={\left({A}^{\text{T}}\right)}^{N-1}\frac{\partial F}{\partial x}\left({x}_{N}^{*}\right)$ (223)

Now consider the system

${y}_{k+1}={D}^{-1}AD{y}_{k}+{D}^{-1}g+{D}^{-1}{e}_{3}{u}_{k}$ (224)

which is conjugate to the ${x}_{k}$ system. Now observe that with

$\stackrel{˜}{H}\left({y}_{k},{v}_{k},{\zeta }_{k}\right)={\zeta }_{k}^{\text{T}}\left(\Lambda {y}_{k}+{D}^{-1}g+{D}^{-1}{e}_{3}{v}_{k}\right)$ (225)

and

$\Lambda ={D}^{-1}AD$ (226)

we have that

${\zeta }_{k}={\left({\Lambda }^{\text{T}}\right)}^{N-k-1}{\zeta }_{N-1}$ (227)

and

${\zeta }_{0}={\left({\Lambda }^{\text{T}}\right)}^{N-1}{\zeta }_{N-1}$ (228)

$\stackrel{˜}{H}\left({y}_{0}^{*},{v}_{0}^{*},{\zeta }_{k}\right)\le \stackrel{˜}{H}\left({y}_{0}^{*},{v}_{0},{\zeta }_{0}\right)$ (229)

is equivalent to

${\zeta }_{0}^{\text{T}}{D}^{-1}B{v}_{0}^{*}\le {\zeta }_{0}^{\text{T}}{D}^{-1}B{v}_{0}$ (230)

So if

$\Delta ={\zeta }_{0}^{\text{T}}{D}^{-1}B<0$ (231)

then maximal chemo therapy is optimal, by a computation like (157) to (163). Here

$\begin{array}{l}{\left({\zeta }_{N-1}\right)}^{\text{T}}=\frac{\partial F}{\partial x}{\left({x}_{N}^{*}\right)}^{\text{T}}\\ =\left(\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right),\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right),-\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\right)\end{array}$ (232)

But this amounts to

$\Delta =\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right){\lambda }_{+}^{N-1}\frac{{D}_{31}}{\text{det}\left(D\right)}$ (233)

$+\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right){\lambda }_{-}^{N-1}\frac{{D}_{32}}{\text{det}\left(D\right)}$ (234)

$\text{ }-\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right){\stackrel{˜}{\lambda }}^{N-1}\frac{{D}_{33}}{\text{det}\left(D\right)}$ (235)

where ${D}_{ij}$ are complements in D. Hence

${D}_{31}=\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\delta \beta -\delta \beta \left(1+{\mu }_{I}-{\lambda }_{-}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)$ (236)

${D}_{32}=-\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\delta \beta +\delta \beta \left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)$ (237)

$\begin{array}{c}{D}_{33}=\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(-\delta \right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}+\delta \left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\end{array}$ (238)

Inserted into (233) to (235) we obtain, that $\Delta$ becomes

$\begin{array}{l}\frac{1}{\text{det}\left(D\right)}\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right){\lambda }_{+}^{N-1}\left(\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\delta \beta \\ -\delta \beta \left(1+{\mu }_{I}-{\lambda }_{-}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\right)\end{array}$ (239)

$\begin{array}{l}+\frac{1}{\mathrm{det}\left(D\right)}\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right){\lambda }_{-}^{N-1}\left(-\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\delta \beta \\ +\delta \beta \left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\right)\end{array}$ (240)

$+\frac{1}{\text{det}\left(D\right)}\left(-\beta \right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right){\stackrel{˜}{\lambda }}^{N-1}\left(\left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(-\delta \right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)$ (241)

$+\delta \left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\right)$ (242)

Now observe that when ${\mu }_{F}={\mu }_{I}=\mu$

$\Delta =\frac{1}{\mathrm{det}\left(D\right)}\left(1+\mu -{\lambda }_{+}\right)\left(1+\mu -{\lambda }_{+}\right)\left(1+\mu -{\lambda }_{-}\right)\left(1+\mu -{\lambda }_{-}\right)\beta \delta \left({\lambda }_{+}^{N-1}-{\lambda }_{-}^{N-1}\right)$ (243)

$=\beta \frac{{\left(\alpha \delta +\beta \sigma \right)}^{2}}{\left({\lambda }_{+}-{\lambda }_{-}\right){\left(\alpha \delta +\beta \sigma \right)}^{2}}\left({\lambda }_{+}^{N-1}-{\lambda }_{-}^{N-1}\right)$ (244)

So maximal chemo therapy is optimal in this case as we have seen above. Now compute

$\Delta =\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\left({\lambda }_{+}^{N-1}-{\lambda }_{-}^{N-1}\right)$ (245)

$+\frac{1}{\text{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(-\left(1+{\mu }_{I}-{\lambda }_{-}\right)\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\left({\lambda }_{+}^{N-1}-{\stackrel{˜}{\lambda }}^{n-1}\right)$ (246)

$+\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\left({\lambda }_{-}^{N-1}-{\stackrel{˜}{\lambda }}^{N-1}\right)$ (247)

Notice that for $k=0,\cdots ,N-2$

$\begin{array}{l}{c}_{k}=\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\cdot \left({\lambda }_{+}^{N-k-1}-{\lambda }_{-}^{N-k-1}\right)<0\end{array}$

when $\mu ={\mu }_{F}={\mu }_{I}$, due to the assumptions. Observe also that (246) = 0, (247) = 0, when $\mu ={\mu }_{F}={\mu }_{I}$, because then $\stackrel{˜}{\lambda }=1+\mu$. Now take $\rho >0$ such that

$\begin{array}{l}|\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(-\left(1+{\mu }_{I}-{\lambda }_{-}\right)\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\\ \begin{array}{c}\text{ }\\ \text{ }\end{array}\cdot \left({\lambda }_{+}^{N-k-1}-{\stackrel{˜}{\lambda }}^{N-k-1}\right)|<-\frac{{c}_{k}}{3}\end{array}$ (248)

and

$\begin{array}{l}|\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\\ \begin{array}{c}\text{ }\\ \text{ }\end{array}\cdot \left({\lambda }_{-}^{N-k-1}-{\stackrel{˜}{\lambda }}^{N-k-1}\right)|<-\frac{{c}_{k}}{3}\end{array}$ (249)

And finally

$\begin{array}{l}-\frac{1}{\mathrm{det}\left(D\right)}\delta \beta \left(1+{\mu }_{F}-{\lambda }_{+}\right)\left(1+{\mu }_{I}-{\lambda }_{+}\right)\left(1+{\mu }_{F}-{\lambda }_{-}\right)\left(1+{\mu }_{I}-{\lambda }_{-}\right)\\ \left({\lambda }_{+}^{N-k-1}-{\lambda }_{-}^{N-k-1}\right)>-\frac{2}{3}{c}_{k}\end{array}$ (250)

But then $\Delta <0$, and maximal chemo therapy is optimal. We have used, that the roots of a polynomial depend continuously on the coefficients of the polynomial, see . For $k=N-1$ notice that $\left(1,0,0\right)\cdot B=0$. So ${u}_{N-1}={g}_{I}^{0}$ is optimal.

When (ii) define the following U below by computing

$\left(\begin{array}{c}\left(1+{\mu }_{F}-a-ib\right)\left(1+{\mu }_{I}-a-ib\right)\\ -\delta \left(1+{\mu }_{I}-a-ib\right)\\ -\sigma \left(1+{\mu }_{F}-a-ib\right)\end{array}\right)$ (251)

$=\left(\begin{array}{c}\left(1+{\mu }_{F}-a\right)\left(1+{\mu }_{I}-a\right)-{b}^{2}+i\left(-\left(1+{\mu }_{I}-a\right)b-\left(1+{\mu }_{F}-a\right)b\right)\\ -\delta \left(1+{\mu }_{I}-a\right)+i\delta b\\ -\sigma \left(1+{\mu }_{F}-a\right)+i\sigma b\end{array}\right)$ (252)

$={v}_{+}={v}_{1}+i{v}_{2}$ (253)

So define U to be

$\left(\begin{array}{ccc}\left(1+{\mu }_{F}-a\right)\left(1+{\mu }_{I}-a\right)-{b}^{2}& -\left(1+{\mu }_{I}-a\right)b-\left(1+{\mu }_{F}-a\right)b& -\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\\ -\delta \left(1+{\mu }_{I}-a\right)& \delta b& \delta \beta \\ -\sigma \left(1+{\mu }_{F}-a\right)& \sigma b& -\alpha \delta +\left(1+\gamma -\stackrel{˜}{\lambda }\right)\left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right)\end{array}\right)$ (254)

Then

${U}^{-1}A\text{\hspace{0.17em}}U=\left(\begin{array}{ccc}a& b& 0\\ -b& a& 0\\ 0& 0& \stackrel{˜}{\lambda }\end{array}\right)$ (255)

Define

$\begin{array}{c}\stackrel{˜}{F}\left(x\right)=\left(\left(1+{\mu }_{F}-a\right)\left(1+{\mu }_{I}-a\right)-{b}^{2}\right){x}_{1}+\left(-b\left(a-1-{\mu }_{I}\right)-b\left(a-1-{\mu }_{F}\right)\right){x}_{2}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}-\beta \left(1+{\mu }_{F}-\stackrel{˜}{\lambda }\right){x}_{3}\end{array}$ (256)

Theorem 3 Suppose (ii) i.e. eigenvalues of A are $a+ib,b\ne 0$. If $N-k-1=2p+1,k=0,\cdots ,N-1$ and

$〈{e}_{1},{A}^{N-k-1}B〉={\frac{\partial \stackrel{˜}{F}}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{A}_{p}& {B}_{p}& 0\\ -{B}_{p}& {A}_{p}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{2p+1}\end{array}\right){U}^{-1}B<0$ (257)

let ${v}_{k}^{*}={g}_{I}^{0}$ and if the reverse inequality holds let ${v}_{k}^{*}=0$. If $N-k-1=2p,k=0,\cdots ,N-1$ and

$〈{e}_{1},{A}^{N-k-1}B〉={\frac{\partial \stackrel{˜}{F}}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{C}_{p}& {D}_{p}& 0\\ -{D}_{p}& {C}_{p}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{2p}\end{array}\right){U}^{-1}B<0$ (258)

let ${v}_{k}^{*}={g}_{I}^{0}$ and if the reverse inequality holds ${v}_{k}^{*}=0$. Then $\left({x}_{k}^{*},{v}_{k}^{*}\right)$ is optimal.

Proof. Define the Hamiltonian

$\stackrel{˜}{H}\left({y}_{k},{v}_{k},{\zeta }_{k}\right)={\zeta }_{k}^{\text{T}}\left({U}^{-1}AU{y}_{k}+{U}^{-1}g+{U}^{-1}B{v}_{k}\right)$ (259)

Then with

$\stackrel{˜}{B}={U}^{-1}AU=\left(\begin{array}{ccc}a& b& 0\\ -b& a& 0\\ 0& 0& \stackrel{˜}{\lambda }\end{array}\right)$ (260)

we find

${\zeta }_{k}={\left({\stackrel{˜}{B}}^{\text{T}}\right)}^{N-k-1}{\zeta }_{N-1}={\left({\stackrel{˜}{B}}^{\text{T}}\right)}^{N-k-1}\frac{\partial \stackrel{˜}{F}}{\partial x}$ (261)

So

$\frac{\partial \stackrel{˜}{H}}{\partial {v}_{k}}={\frac{\partial \stackrel{˜}{F}}{\partial x}}^{\text{T}}{\stackrel{˜}{B}}^{N-k-1}{U}^{-1}B={\frac{\partial \stackrel{˜}{F}}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{A}_{p}& {B}_{p}& 0\\ -{B}_{p}& {A}_{p}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{2p+1}\end{array}\right){U}^{-1}B$ (262)

when $N-k-1=2p+1,k=0,\cdots ,N-1$ and when $N-k-1=2p,k=0,\cdots ,N-1$

${\frac{\partial \stackrel{˜}{F}}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{C}_{p}& {D}_{p}& 0\\ -{D}_{p}& {C}_{p}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{2p}\end{array}\right){U}^{-1}B$ (263)

Optimality follows from a computation like (157) to (163).

Theorem 4 Suppose (i) and the eigenvaues of A are real and distinct. If

$〈{e}_{1},{A}^{N-k-1}B〉={\frac{\partial F}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{\lambda }_{+}^{N-k-1}& 0& 0\\ 0& {\lambda }_{-}^{N-k-1}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{N-k-1}\end{array}\right){D}^{-1}B<0$ (264)

let ${v}_{k}^{*}={g}_{I}^{0}$ and if the reverse inequality holds let ${v}_{k}^{*}=0$. Then $\left({x}_{k}^{*},{v}_{k}^{*}\right)$ is optimal.

Proof. Define the Hamiltonian

$\stackrel{˜}{H}\left({y}_{k},{v}_{k},{\zeta }_{k}\right)={\zeta }_{k}^{\text{T}}\left({D}^{-1}AD{y}_{k}+{D}^{-1}g+{D}^{-1}B{v}_{k}\right)$ (265)

Then with

$\Lambda ={D}^{-1}AD$ (266)

we get

${\zeta }_{k}={\left({\Lambda }^{\text{T}}\right)}^{N-k-1}{\zeta }_{N-1}={\left({\Lambda }^{\text{T}}\right)}^{N-k-1}\frac{\partial F}{\partial x}$ (267)

Then the partial derivative of $\stackrel{˜}{H}$ with respect to ${v}_{k}$ is

${\frac{\partial F}{\partial x}}^{\text{T}}\left(\begin{array}{ccc}{\lambda }_{+}^{N-k-1}& 0& 0\\ 0& {\lambda }_{-}^{N-k-1}& 0\\ 0& 0& {\stackrel{˜}{\lambda }}^{N-k-1}\end{array}\right){D}^{-1}B$ (268)

,

Optimality follows from a computation like (157) to (163).

6. Summary

In the present paper we applied the discrete Pontryagin Minimal Principle to a discrete model T of cancer growth and the Pontryagin Minimal Principle to an affine vector field that generates T. When $\mu ={\mu }_{F}={\mu }_{I}$ and the eigenvalues of T are real and distinct, maximal chemo therapy is optimal for the discrete model, while this is not necessarily so when the eigenvalues of A are $1+\mu ,a+ib,b\ne 0,1+\mu >0$.

For the affine vector field that generates T, we have proven similar statements, when $\mu ={\mu }_{F}={\mu }_{I}$. Maximal chemo therapy is optimal, when the eigenvalues of A are real, positive and distinct and this is not necessarily so, when there are imaginary eigenvalues. We finally considered what happens in the discrete model, when ${\mu }_{F}\ne {\mu }_{I}$. In particular we have derived an optimal strategy to give chemo or immune therapy.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

Cite this paper

Larsen, J.C. (2019) Optimal Control of Cancer Growth. Applied Mathematics, 10, 173-195. https://doi.org/10.4236/am.2019.104014

References

1. 1. Chr, J. (2016) Larsen Models of Cancer Growth. Journal of Applied Mathematics and Computing, 53, 615-643.

2. 2. Seierstad, A. and Sydsaeter, K. (1988) Optimal Control Theory with Economic Applications North Holland.

3. 3. Goodwin, G. (2005) Constrained Control and Estimation. An Optimization Approach. Springer Verlag, Berlin.

4. 4. Swierniak, A., Kimmel, M. and Smieja, J. (2009) Mathematical Modelling as a Tool for Planning Anti Cancer Therapy. European Journal of Pharmacology, 625, 108-121. https://doi.org/10.1016/j.ejphar.2009.08.041

5. 5. Laird, A.K. (1964) Dynamics of Cancer Growth. British Journal of Cancer, 18, 490-502. https://doi.org/10.1038/bjc.1964.55

6. 6. Adam, J.A. and Bellomo, N. (1997) A Survey of Models for Tumor-Induced Immune System Dynamics. Birkhauser, Boston.

7. 7. Geha, R. and Notarangelo, L. (2012) Case Studies in Immunology: A Clinical Companion. Garland Science, Hamden.

8. 8. Marks, F., Klingmüller, U. and Müller-Decker, K. (2009) Cellular Signal Processing. Garland Science, Hamden.

9. 9. Molina-Paris, C. and Lythe, G. (2011) Mathematical Models and Immune Cell Biology. Springer Verlag, Berlin. https://doi.org/10.1007/978-1-4419-7725-0

10. 10. Murphy, K. (2012) Immunobiology. 8th Edition, Garland Science, Hamden.

11. 11. Rees, R.C. (2014) Tumor Immunology and Immunotherapy. Oxford University Press, Oxford.

12. 12. Larsen, J.C. (2016) Hopf Bifurcations in Cancer Models. JP Journal of Applied Mathematics, 14, 1-31.

13. 13. Larsen, J.C. (2017) A Mathematical Model of Adoptive T Cell Therapy. JP Journal of Applied Mathematics, 15, 1-33.

14. 14. Larsen, J.C. (2016) Fundamental Concepts in Dynamics.

15. 15. Larsen, J.C. (2017) The Bistability Theorem in a Cancer Model. International Journal of Biomathematics, 11, Article ID: 1850004.

16. 16. Larsen, J.C. (2016) The Bistability Theorem in a Model of Metastatic Cancer. Applied Mathematics, 7, 1183-1206. https://doi.org/10.4236/am.2016.710105

17. 17. Larsen, J.C. (2017) A Study on Multipeutics. Applied Mathematics, 8, 746-773. https://doi.org/10.4236/am.2017.85059

18. 18. Larsen, J.C. (2017) A Mathematical Model of Immunity. JP Journal of Applied Mathematics.

19. 19. Larsen, J.C. (2018) Models of Cancer Growth Revisited. Applied Mathematics, 9, Article ID: 84308.

20. 20. Uspensky (1948) Theory of Equations. McGraw-Hill, New York.

21. 21. Alexendarian, A. (2013) On Continuous Dependence of Roots of Polynomials on Coefficients.