% SIMIODE-TeacherStudent-Template.tex and simiode.cls for LaTeX2e
% Version 22 November 2016.
%
% This template is to be used with the  class file "simiode.cls".
%
% The latest version of this class and template files
% can be found at
%
%       <http://www.SIMIODE.org> in the Author section.
%
% A "Quick Start Guide" appears below.
%
% This file must be processed with version LaTeX2e or higher.
% Appreciation to Mark Frantz, IUPUI; John Thoo, Yuba College; and Eric Sullivan, Carroll College.
%
\documentclass{simiode}
%
% With this version of the tex file and the associated SIMIODE.cls file you can produce BOTH a Student Version and a Teacher Version of your modeling scenario or technique narrative.
%
% Student Version contains the title STUDENT VERSION along with the NAME OF SCENARIO, with no abstract, no keywords, no tags, no author identifiers, and no COMMENTS section. The Student Version contains only the STATEMENT of the modeling scenario only with appropriate figures, data, instructions, and perhaps a bibliography appropriate for student resourcing, but nothing in the bibliography which would "give away" the nature of the modeling scenario, etc.
%
%The Teacher Version contains the material in the Student Version, namely the title NAME OF SCENARIO, now headed with   TEACHER VERSION, but also full author  identification, abstract, key words, tag words, and, following the STATEMENT, a section of COMMENTS which can be discussion of presentation  ideas, "solutions," techniques, etc.
%
% Authors can prepare a full Teacher Version using the notions immediately below  to create a pdf for review, but in review please eliminate all information about the author as we use a double-blind referee system. Such information will be restored when accepted and published in SIMIODE at www.simiode.org.
%
%***********   Toggle the following line to turn the teacher version on (1) or off (0). ***********
%
\def\TeacherEdition{0} % Toggle this line's value.
\newcommand{\teacher}[1]{\ifnum\TeacherEdition=1 #1 \fi}
%
% The \teacher{ } command wraps any text that will appear only in the teacher version of
% the document. See this file for examples. This has been done for you in the document below.
% Rename the pdf after you compile each final version.
%
\usepackage{amsmath,amssymb}
\usepackage{color}
\usepackage{graphicx}
\usepackage{hyperref}
\usepackage{url}
%
\renewcommand{\baselinestretch}{1.3}
%
% NOTE regarding macros %%
%
% Any macros defined for your paper should be contained in
% the top matter. Likewise, any environment definitions such
% as \newtheorem or \newenvironment should also go in the
% top matter.
%
% Do not \input your macros as separate files within your
% final version because more files make it harder for the
% editors and users to keep track of and/or modify and
% append your paper.
%
% Top matter is above and the beginning of the document below.
%
\begin{document}

\newpage
%
% You only need to change the NAME OF SCENARIO in the next line and delete the footnote.
% See above for instructions on toggling from Teacher Version to Student Version.
%
\title{\ifnum\TeacherEdition=1 TEACHER VERSION \else STUDENT VERSION \fi\\Modeling Cancer Growth with Differential Equations}



\markboth{Modeling Cancer Growth}{Modeling Cancer Growth}

\author{Jue Wang\\
Union College\\
Schenectady NY 12308 USA\\
wangj@union.edu
}

\makeStitlePDFLaTex


\begin{abstract}
This module guides students in the use of differential equation models to predict cancer growth and optimize treatment outcomes. Several classical models for cancer growth are studied, including  exponential, power law, Bertalanffy, logistic, and Gompertz. Students solve first-order differential equations using separation of variables and/or partial fractions. They examine the behaviors of the equations and solutions through qualitative techniques. Furthermore, students evaluate how cancer treatments affect tumor growth. Real cancer data are provided to give options for data fitting and prediction of cancer growth.
%It addresses some main topics in differential equations, such as initial value problems, analytic solutions, qualitative approaches, and numerical approximation. The techniques include solving first-order differential equations using separation of variables and partial fractions, equilibria and phase line, slope fields and solution graphs, stability analysis, and Euler method. 

Through step-by-step investigation, conjecturing, predicting, and analyzing, students discover how their knowledge can be used to address complex and real problems, and improve problem solving ability and mathematical reasoning skills.

\end{abstract}

\section*{SCENARIO DESCRIPTION}
Cancer is approaching  heart disease as the leading cause of death in the United States \cite{ASCO2017}. Knowledge of tumor growth is important in cancer screening planning and treatment. Effective models of tumor growth not only help researchers understand the dynamics of cancer growth and evaluate screening strategies, but they also enable clinicians to predict the efficacy of cancer therapies and optimize treatment protocols.

In this project, you will investigate several classical models for cancer growth and examine the behaviors of these equations and their solutions through qualitative techniques. Furthermore, you will evaluate how the predictions of tumor growth change when chemotherapy or other treatments are introduced, such as, how does the maximum tumor size change? What is the minimum effective concentration of the drug needed? Analysis and data fitting may be carried out to assess the models' ability to predict future growth of tumors and shed some light on the development of resistance to cancer treatment.

\section{Modeling scenarios}
The growth of cancer can be modeled through first-order ordinary differential equations by examining the change of tumor volume $V$ over time $t$, given an initial volume $V(0)=V_0$. This initial condition is used in the following models. The parameters in the models are positive.

\subsection{Exponential growth}
A simple exponential model describes the early stages of tumor growth, where the growth is proportional to the tumor population. Cells divide and create two daughter cells at a constant rate that is independent of tumor size. The model is given by (1). This model was first applied to cancer growth in 1956 \cite{Collins1956}.
\begin{equation}
\frac{dV}{dt}=aV
\end{equation}

The tumor doubling time is introduced to measure how fast a tumor grows and quantify the growth rate, $$DT=\frac{\ln 2}{a}.$$

The exponential model predicts early growth well. However, with depletion of nutrients and cell death, exponential growth is not valid for the long-term growth of solid tumors. Thus we investigate several alternative models below.\\ \\
{\bf Exercise 1.} Solve (1) for tumor volume $V(t)$. Determine the maximum tumor size.

\subsection{Power law}
A generalization of the exponential model is described by the power-law equation (2), introduced by Mendelsohn in 1963 \cite{Mendelsohn1963}.
\begin{equation}
\frac{dV}{dt}=aV^b
\end{equation}
\\
{\bf Exercise 2.} Solve (2) for tumor volume $V(t)$ for $b<1$. Determine the maximum tumor size.

\subsection{Von Bertalanffy}
The Bertalanffy model balances synthesis and destruction \cite{Bertalanffy1957}. Assume that growth is proportional to surface area, since nutrient enters through the surface, and that death is proportional the tumor size. The model is given by (3), also known as the surface rule model.
\begin{equation}
\frac{dV}{dt}=aV^{2/3}-bV
\end{equation}
This model has been successfully applied to describe human tumor growth. \\ \\
{\bf Exercise 3.} The exact solution of  (3) is given below. Determine the maximum tumor size.
$$V(t)=\left[\frac{a}{b}+\left(V_0^{1/3}-\frac{a}{b}\right)e^{-bt/3}\right]^3.$$

\subsection{Logistic}
The logistic model of population growth was first published by Pierre Verhulst in 1845 \cite{Verhulst1845}. The model assumes a linear decrease of the relative growth rate with population size. The maximum size is limited by a carrying capacity $K$.
\begin{equation}
\frac{dV}{dt}=aV\left(1-\frac{V}{K}\right)
\end{equation}
\\
{\bf Exercise 4.} Solve (4) for tumor volume $V(t)$ for $a=1$ and $K=1$. Determine the maximum tumor size.

\subsection{Gompertz}
The Gompertz model was created by Benjamin Gompertz in 1825 to explain human mortality curves and determine the value of life insurances \cite{Gompertz1825}. It was applied to organism growth a hundred years later. The Gompertz model exhibits an exponential decay of the growth rate. It has been successfully used to model breast and lung cancer growth.
\begin{equation}
\frac{dV}{dt}=r(t)V(t), \quad \frac{dr}{dt}=-br(t)
\end{equation}

From the above equations, we have
$$
\frac{d}{dt}\ln V=\frac{1}{V}\frac{dV}{dt}=r(t)=-\frac{1}{b}\frac{dr}{dt}.
$$

Thus,
$$
\ln V=\frac{-r(t)+a}{b} \quad \textnormal{for some constant } a.
$$

Solving for $r(t)$ leads to $r(t)=a-b\ln V$. Substituting this into (5) we obtain an alternative and more popular form of the Gompertz model.
\begin{equation}
\frac{dV}{dt}=V(a-b\ln V)
\end{equation}
\\
{\bf Exercise 5.} Rewrite (5) as $dV/dt=k e^{-bt}V$. Solve it for tumor volume $V(t)$. Determine the maximum tumor size.

\subsection{Modeling summary}
In the case of multicellular tumor spheroids, cancer growth follows a sigmoid curve with three phases: initial exponential growth, linear growth, and plateau. The Bertalanffy, logistic, and Gompertz models reflect the sigmoid nature of growth. They are summarized in Figure 1, originating from the generalized two-parameter model of growth for biological systems \cite{Marusic1996}. Note that these models are autonomous differential equations.

\begin{figure}[htb]
\includegraphics[width=3.7in]{growthmodel}
\centering
\caption{Nesting scheme for empirical models of growth \cite{Marusic1996}.}
\end{figure}

\section{Model fitting and analysis}
In the second part of this activity, you will fit  models to real cancer data to predict the best fit model for tumor growth. The fitting is performed by minimizing the sum of squared residuals (SSR) or normalized mean square error (NMSE). The NMSE is independent of measurement units and scale. A value of NMSE close to zero indicates a good fit of the model.
\begin{equation}
SSR=\sum_i (y_i-\hat{y}_i)^2
\end{equation}
\begin{equation}
NMSE=\frac{\sum_i (y_i-\hat{y}_i)^2}{\sum_i y_i^2}
\end{equation}

The $\/{y_i\/}$'s are tumor volume measurements at time points $t_i$; the $\/{\hat{y}_i\/}$'s are tumor volume estimates predicted by the model at the same time points.

\subsection{Fibroblast tumor cell}
The dataset provided in Table 1 is taken from the Chinese hamster V79 fibroblast tumor cell line \cite{Marusic1994}. It consists of 45 measurements of volumes ($10^9 \mu m^3$) during the time period of 60 days.

\begin{table}[htb]
\begin{center}

\begin{tabular}{ | c | c || c | c || c | c || c | c || c | c | }
\hline
	{\bf t} (days) & \bf V ($10^9 \mu m^3$) & \bf t & \bf V & \bf t & \bf V & \bf t & \bf V & \bf t & \bf V \\ \hline
	3.46 & 0.0158 & 12.39 & 0.4977 & 24.33 & 3.2046 & 35.20 & 5.9668 & 48.29 & 7.0694 \\ \hline
	4.58 & 0.0264 & 13.42 & 0.6033 & 25.58 & 4.5241 & 36.34 & 6.6945 & 49.24 & 7.4971 \\ \hline
	5.67 & 0.0326 & 15.19 & 0.8441 & 26.43 & 4.3459 & 37.29 & 6.6395 & 50.19 & 6.9974 \\ \hline
	6.64 & 0.0445 & 16.24 & 1.2163 & 27.44 & 5.1374 & 38.50 & 6.8971 & 51.14 & 6.7219 \\ \hline
	7.63 & 0.0646 & 17.23 & 1.4470 & 28.43 & 5.5376 & 39.67 & 7.2966 & 52.10 & 7.0523 \\ \hline
	8.41 & 0.0933 & 18.18 & 2.3298 & 30.49 & 4.8946 & 41.37 & 7.2268 & 54.00 & 7.1095 \\ \hline
	9.32 & 0.1454 & 19.29 & 2.5342 & 31.34 & 5.0660 & 42.58 & 6.8815 & 56.33 & 7.0694 \\ \hline
	10.27 & 0.2183 & 21.23 & 3.0064 & 32.34 & 6.1494 & 45.39 & 8.0993 & 57.33 & 8.0562 \\ \hline
	11.19 & 0.2842 & 21.99 & 3.4044 & 33.00 & 6.8548 & 46.38 & 7.2112 & 59.38 & 7.2268 \\ \hline
\end{tabular}
\caption{Chinese hamster V79 fibroblast tumor}

\end{center}
\end{table}

{\bf Evaluate the models and perform the following tasks.}
\begin{enumerate}
\item {\bf The best model for cancer growth.} \\
Consider the three models $V'=f(V)$: {\bf Bertalanffy}, {\bf logistic}, and {\bf Gompertz}. Use model equations given in Figure 1. Note that $K=a/b$ for logistic model. Find the optimal values of parameters $a$ and $b$ by minimizing the NMSE. Which model best describes the cancer growth? You may use Excel ``Solver" function or Matlab or any other software of your choice. \\

Hint: In Excel, the model estimate can be evaluated using Euler's method to numerically solve the differential equations:
\begin{equation}
V_i=V_{i-1}+f(V_{i-1})\cdot(t_i-t_{i-1}), \quad i=1,2,\ldots, 45\, .
\end{equation}

\noindent A helpful guide to Excel ``Solver" function is provided in \cite{ExcelSolver}.\\

\item The {\bf optimal parameter values} and NMSE obtained from above are:
\begin{center}
\begin{tabular}{ll}
{\bf Bertalanffy} & $a=0.4340, \quad b=0.2158, \quad NMSE=0.0089$ \\
{\bf Logistic} & $a=0.3389, \quad b=0.0489, \quad NMSE=0.0138$ \\
{\bf Gompertz} & $a=0.2375, \quad b=0.1179, \quad NMSE=0.0049$
\end{tabular}
\end{center}
The Gompertz model gives the best fit to the tumor growth data with an NMSE of 0.49\%.

\begin{enumerate}
\item Find the equilibria of the three models using these values. Plot the three growth rate functions $f(V)$ in one graph. Label graph including proper units.
\item Sketch the phase line for each model. Identify the types of equilibria (sink, source, node) and their stabilities.

\item Hand sketch the slope fields and solution graphs in separate planes.
\item In Excel, or Matlab, or other softwares, plot the data measurements given in Table 1 and three model fitting curves (9) in one graph. Label your graph including proper units.
\item Describe the behavior of the solutions.
\end{enumerate}
\smallskip
\item {\bf Prediction of cancer growth.} \\
Repeat the task in Problem 1, but only use the first third of the data to find the optimal values of $a$ and $b$. How well do the models predict cancer growth? Which model best predicts cancer growth? Assess the difference in model predictions. You may compute the NMSE and plot the solution graphs. Can the model prediction be improved with more data (such as half of the data)? How much improvement can we see? \\

\item {\bf With cancer treatments} (e.g., chemotherapy, immunotherapy, or drug treatment).\\
The prediction will be affected when cancer treatment is added to the model. It is important to examine this effect because the cancer growth model is often used as a basis for assessing the efficacy of cancer treatments. We model a simple implementation of a cancer treatment by subtracting $cV$ from each model equation for some constant drug supply $c$ \cite{Murphy2016}.
\begin{enumerate}
\item What is the maximum tumor volume for all five models with cancer treatment? Do not plug in the parameter values. Express your answers in terms of $a, b, c$. (Hint: You don't need to solve the differential  equations. You may find the equilibria and use stability analysis.)
\item What is the minimum amount of drug needed to eradicate the tumor in each model? (i.e. to achieve a negative initial growth rate of the tumor.) Express your answers in terms of $a, b, c, V_0$. Can the tumor be eradicated?
\end{enumerate}

\end{enumerate}

\subsection{Mouse skin carcinogenesis}
The second dataset provided in Table 2 is taken from mice tumors that arise in situ from normal skin as a result of chemical mutagenesis, and develop naturally in the context of a full immune system, thus more closely resembling the complexity and dynamics of tumorigenesis {\it in vivo} \cite{Loizides2015}. It consists of 24 measurements during the time period of 11 weeks.

\begin{table}[htb]
\begin{center}

\begin{tabular}{ | c | c || c | c || c | c || c | c | }
\hline
	{\bf t} (days) & \bf V ($mm^3$) & {\bf t} (days) & \bf V ($mm^3$) & {\bf t} (days) & \bf V ($mm^3$) & {\bf t} (days) & \bf V ($mm^3$) \\ \hline
	0    & 1.0410  & 21  &  35.1929   & 42  &  87.7625    & 62	& 407.5250 \\ \hline
	3    & 4.3156  & 24  &  32.3758   & 45  &  98.7135    & 65	& 303.3899 \\ \hline
	7    & 2.8523  & 28  &  36.9683   & 49  &  131.8011  & 69	& 393.7301 \\ \hline
	10  & 4.1246  & 31  &  62.4162   & 52  &  147.4074  & 72	& 337.5277 \\ \hline
	14  & 5.9634  & 35  & 135.3351  & 56  &  267.1758  & 76	& 376.4903 \\ \hline
	17  & 8.2353  & 38  & 137.0422	  & 59  &  421.1538  & 78	& 395.8532 \\ \hline
\end{tabular}
\caption{Mouse tumor CM.37 T.1}
\end{center}
\end{table}

Perform the same tasks as in problems 1 \& 2 in Sec. 2.1. The optimal fit values and NMSE are found below. Which model gives the best fit? Is it the same or different? If different, what are possible reasons? Examine and compare the two dataset plots. Has the second dataset reached the third phase, {\it i.e.}, plateau?
\begin{center}
\begin{tabular}{lll}
{\bf Bertalanffy} & $a=0.4947, \quad b=0.0465$, & $NMSE=0.0514$ \\
{\bf Logistic} & $a=0.1402, \quad b=0.000357$, & $NMSE=0.0388$ \\
{\bf Gompertz} & $a=0.2619, \quad b=0.0413$, & $NMSE=0.0463$
\end{tabular}
\end{center}

\section{Conclusions}
Write a paragraph summarizing your analysis, findings, and recommendations.


%\section*{Submission Version}

%For submission to the manuscript management system, \url{https://simiode.expressacademic.org} authors should submit a Teacher Version with all author identifying information stripped out. This should be a pdf version of the material.

%Authors may also submit supporting documents such as pdf's of code files, e.g., MatLab, Mathematica, Maple, SAGE, etc. and spreadsheets and data files.


%Below is a set of References which appears in the Student Version. Note, there is a second set of References one can use for the Teacher Version.

\begin{thebibliography}{98}

\bibitem{ASCO2017} ASCO, The state of cancer care in America, 2017: A report by the American Society of Clinical Oncology. 2017. {\it J. Oncol. Pract.} 13(4): e353--e394.

\bibitem{Collins1956} Collins, V.P., R.K. Loeffler and H. Tivey. 1956. Observations on growth rates of human tumors. {\it Am. J. Roentgeol. Radium Ther. Nucl. Med.\/} 76(5): 988--1000.


\bibitem{Bertalanffy1957} Bertalanffy, L. von. 1957. Quantitative laws in metabolism and growth. {\it Q. Rev. Bioi.\/} 32: 217--231.

\bibitem{ExcelSolver} Nonlinear model fitting using Excel. \href{https://youtu.be/3Fd4ukzkxps}{https://youtu.be/3Fd4ukzkxps}. Accessed 29 July 2018.




\bibitem{Gompertz1825} Gompertz, B. 1825. On the nature of the function expressive of the law of human mortality, and on a new method of determining the value of life contingencies. {\it Phil. Trans. Roy. Soc.\/}, 27(5):  13--85.


\bibitem{Loizides2015}Loizides,  C., D. Iacovides, M.M. Hadjiandreou, G. Rizki, A. Achilleos, K. Strati and G.D. Mitsis. 2015.  Model-based tumor growth dynamics and therapy response in a mouse model of de novo carcinogenesis. {\it PLOS ONE\/}. 10(12): 1--18.

\bibitem{Marusic1996} Maru\u{s}i\'{c}, M. 1996. Mathematical models of tumor growth, {\it
Lect. Math. Colloquium Croatian Math. Soc.\/}, Osijek. 175--192, June 7.

\bibitem{Marusic1994} Maru\u{s}i\'{c}, M., \u{Z}. Bajzer, J.P. Freyer and S.Vuk-Pavlovi\'{c}. 1994.  Analysis of growth of multicellular tumour spheroids by mathematical models. {\it
Cell Prolif.\/} 27: 74--94.

\bibitem{Mendelsohn1963} Mendelsohn, M. L. 1963. Cell proliferation and tumor growth. {\it Cell Proliferation\/}. Oxford: Blackwell Scientific Publications. 190--210

\bibitem{Murphy2016} Murphy, H., H. Jaafari and H. Dobrovolny,. 2016. Differences in predictions of ODE models of tumor growth: a cautionary example. {\it BMC Cancer\/}, 16: 163--172.

\bibitem{Verhulst1845} Verhulst, P.F.  1845. Recherches math\'{e}matiques sur la loi d'accroissement de la population. {\it Nouv. m\'{e}m. de l'Academie Royale des Sci. et Belles-Lettres de Bruxelles\/}. 18: 1--41.

\end{thebibliography}

\teacher{
\section*{NOTES FOR TEACHERS}

This module is designed for an introductory differential equations course. It can be used early in the course as an individual assignment or group project, or combination (part one individual and part two group work). Students may use Excel, Matlab, or other softwares to perform data fitting and parameter optimization. The instructor has the option to use only the first dataset or both, and/or assign the data fitting part as bonus or optional, since the optimal parameter values and NMSE are provided. A simple example of data fitting may be demonstrated in class, such as in \cite{ExcelSolver}. Students may spend 8-10 hours on the project. 

The Gompertz model gives the best fit to the first dataset, but the logistic model gives the best fit to the second data. The second dataset has not reached the third phase (plateau). The choice of cancer growth model may lead to different prediction outcomes.

More datasets are available upon request.

\subsection*{Partial Solutions:}
1.

{\bf Exercise 1.} $V(t)=V_0 e^{at}, \quad \lim_{t\rightarrow +\infty} V(t) = +\infty$.

{\bf Exercise 2.} $V(t)=(V_0^{1-b}+(1-b)at)^{1/(1-b)}, \quad \lim_{t\rightarrow +\infty} V(t) = +\infty$.

{\bf Exercise 3.} $\lim_{t\rightarrow +\infty} V(t) = (a/b)^3$.

{\bf Exercise 4.} $V(t)=\frac{V_0 K}{V_0+(K-V_0)e^{-at}}$ for general $a$ and $K$, \quad $\lim_{t\rightarrow +\infty} V(t) = K$.

{\bf Exercise 5.} $V(t)=V_0e^{\frac{k}{b}(1-e^{-bt})}$, \quad $\lim_{t\rightarrow +\infty} V(t) = V_0e^{\frac{k}{b}}$.

An alternative approach is to solve for the non-zero equilibrium of (6): $V_\textnormal{max}=e^\frac{a}{b}$.\\ \\
2.1: 2(b)
\begin{center}
\begin{tabular}{lll}
{\bf Bertalanffy} & $V=0$ (node, semi-stable), & $V=(a/b)^3=8.1342$ (sink, stable) \\
{\bf Logistic} & $V=0$ (source, unstable), & $V=a/b=K=6.9305$ (sink, stable) \\
{\bf Gompertz} & & $V=e^\frac{a}{b}=7.4964$ (sink, stable)
\end{tabular}
\end{center}

\newpage
2.1: 2(d)
\begin{figure}[htb]
\includegraphics[width=5.5in]{data1}
\centering
\caption{Dataset in 2.1.}
\end{figure}

2.2
\begin{figure}[htb]
\includegraphics[width=5.5in]{data2}
\centering
\caption{Dataset in 2.2.}
\end{figure}

%\begin{thebibliography}{99}
%\end{thebibliography}
}


\end{document}
