% 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
% 
%       \leq 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\\Enabling Epidemic Exploration Education}




\markboth{Enabling Epidemic Exploration Education}{Enabling Epidemic Exploration Education}

\author{Brian Winkel,
Director, SIMIODE\\ 
Chardon OH 44024 USA\\
BrianWinkel@simiode.org
}

\makeStitlePDFLaTex


\begin{abstract}
We offer several strategies for estimating parameters in models of epidemics, one using a Michaelis-Menten saturation infected rate.
\end{abstract}


\section*{SCENARIO DESCRIPTION}

While working with the principal author of \cite{MillerHelms2010} to prepare the essence \cite{Switkes2003} of a modeling scenario on the basic SIR model for epidemics for wider use than as a project in the author's class we became aware of several interesting possibilities. Since  \cite{MillerHelms2010} presents a modeling opportunity with data obtained and cited in \cite{CDC1978, Murray1993} we explored the several routes to parameter estimation in the SIR model with respect to the data pictured in \cite{CDC1978} and offered in Table~\ref{table:data}.

\begin{table}[!h]
\begin{center}
\vspace*{.2 true cm}
\small
\begin{tabular}{|l|c|c|c|c|c|c|c|c|c|c|c|c|c|c|} \hline
Time (in days)	& 0& 	1& 	2& 	3	& 4& 	5& 	6& 	7& 	8& 	9& 	10	& 11& 	12	& 13\\ \hline
\# Bedridden Boys& 	1	& 3& 	25	& 72	& 222	& 282	& 256	& 233& 	189	& 123& 	70	& 25	& 11	& 4\\ \hline
\end{tabular}


\caption{\label{table:data}Total Number of Bedridden Boys on Day $t$.}
\end{center}
\end{table}

Further, we move to a new SIR model in which we employ a Michaelis-Menten form of a saturated incidence rate between infectives and susceptibles. We suggest several strategies for estimating parameters in this newer model and challenge the reader  to compare our results to the more traditional SIR model efforts.

From \cite{MillerHelms2010} we quote

\begin{quote}
A boarding school is a somewhat closed community; all the students and teachers tend to live on or near campus, and the students do not regularly interact with people not in the boarding school community. Table~\ref{table:data} gives data for a flu outbreak at a boarding school in England. No one died in this flu outbreak as it was not Spanish Flu. These data were compiled from the Communicable Disease Surveillance Center \cite[p.~587]{CDC1978}  and are given as an example by Murray in his book on mathematical biology \cite{Murray1993}. The data values were extracted from the graph by \cite{Switkes2003}.
\end{quote}

\noindent  From \cite{Miller2014}  we quote

\begin{quote}
There are 763 boys in the school. From Table~\ref{table:data}  we can see that the initial number of bedridden students is one. Let us assume that the number of boys who are bedridden
on day $t$ (the data given) is the number of students infected on day $t$, i.e.  no boys are bedridden for any other reason during this period of time and all boys infected will be bedridden. One may wish to test the effects of this last assumption on the model in the analysis phase. Another, perhaps more realistic assumption, is that the bedridden boys represent a fixed percentage of the infected students. This assumption, which excludes teachers and staff from the population, implies that the number of students infected on day zero is one, i.e.~$I(0) = 1$, where $I(t)$ is the number of students infected on day $t$.
\end{quote}

First, we offer a traditional system of differential equations (\ref{eq:SIR1}) for the spread of disease from susceptible, to infected, to recovered individuals. This can be depicted in the diagram in Figure~\ref{fig:6-019-EpidemicArticleArrowSIRDiagram}.

\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=9.5cm]{6-019-EpidemicArticleArrowSIRDiagram}



\parbox{11. true cm}{\caption{\label{fig:6-019-EpidemicArticleArrowSIRDiagram} Picture for a traditional SIR model of an epidemic with connecting rates $a\cdot S(t)\cdot I(t)$ and $b\cdot I(t)$. }}
\vspace*{-.3 true cm}
\end{figure}




Our system of differential equations in (\ref{eq:SIR1})  is referred to as the {\it SIR model for the spread of disease\/} or an epidemic model.
\begin{eqnarray}\label{eq:SIR1}\nonumber
S'(t) &=& —a\cdot S(t)\cdot I(t)\\ 
V(t) &=& a\cdot S(t)\cdot I(t) - b\cdot I(t)\\ \nonumber
R'(t)& =& b\cdot I(t)
\end{eqnarray}

\noindent where $S(t)$ is the number of susceptible students at time $t$, $I(t)$ is the number of infected and infectious students at time $t$, and $R(t)$ is the number of the infected students who have recovered at time $t$.

Here $a$ is the infection rate or proportion with units 1/((infected number)(time)), while $a\cdot S(t)\cdot I(t)$ represents the number of possible interactions between the susceptible students and the infected students which result in an increase in the infected student population per unit time. $b$ is the recovery rate with units (recovered number) /((infected number)(time)).

Second, we offer a model (\ref{eq:SIR2}) for the spread of disease in which there is a form of a saturated incidence rate between infectives and susceptibles. 

\begin{eqnarray}\label{eq:SIR2}\nonumber
S'(t) &=& -a\frac{S(t)}{K_m + S(t)}I(t)\\
I(t)& =& a\frac{S(t)}{K_m + S(t)}I(t) - b\cdot I(t)\\ \nonumber
R'(t) &=& b\cdot I(t)
\end{eqnarray}



This can be depicted in the diagram in Figure~\ref{fig:6-019-EpidemicArticleArrowSIRDiagramMM} 
 where the variables in (\ref{eq:SIR2}) are as in (\ref{eq:SIR1})..


\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=9.5cm]{6-019-EpidemicArticleArrowSIRDiagramMM}



\parbox{11. true cm}{\caption{\label{fig:6-019-EpidemicArticleArrowSIRDiagramMM} Picture for an SIR model of an epidemic which models saturation or satiation of the incidence rate between infectives and susceptibles with connecting rates $a\frac{S(t)}{K_m + S(t)}  I(t)$ and $b\cdot I(t)$.}}
\vspace*{-.3 true cm}
\end{figure}






 

Here $a$ is the infection rate or proportion with units (susceptible number)/((infected number)(time)), while $a \frac{S(t)}{K_m + S(t)}  I(t)$ represents the number of possible interactions between the susceptible students and the infected students which result in an increase in the infected student population per unit time. In this model $K_m$ represents the number of susceptible students at which a given infected student will infect a susceptible student at one half the maximum infection rate, which is $a$. See Figure~\ref{fig:6-019-EpidemicArticleHolling}.


\begin{figure}[!h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=10.5cm]{6-019-EpidemicArticleHolling}



\parbox{12. true cm}{\caption{\label{fig:6-019-EpidemicArticleHolling} The functional response $\frac{S(t)}{K_m + S(t)}S(t) I(t)$ for $K_m = 40$. The vertical line depicts the susceptible population level and the horizontal line points to the 50\% or 0.50 level of infection corresponding to that susceptible population level, exactly one half the maximum infection value of 1 in this case.}}
\vspace*{-.3 true cm}
\end{figure}

\teacher{
\section*{What's in it for student and teacher?}

Before proceeding to the details of the models and analyses let us consider what is in this venture for students and teachers.

Students need to be motivated to learn mathematics and the full modeling process from soup to nuts is one way to do this, a way we have found very successful during our teaching career. Modeling offers students a chance to seek out and learn mathematics, often in a just in time manner, to help advance their cause. Students put this mathematics in a context which is immediate and meaningful for them and they see the payoff in modeling success. They also see failure in both their model and their ability to apply the mathematics they believe is appropriate, e.g., not every problem can be put into a set of linear constraints with a single objective function to be slain by the Simplex Algorithm of linear programming lore! However, in success students see mathematics (and statistics) as part of the solution and validation of other models. Most of the time we do not permit students the ``luxury'' of this entire process. Indeed, all too often we only let our students see the mathematics out of context and without context.

No matter what level of mathematics is under study the learning can always be enhanced by putting that mathematics in context. In \cite[p. 92]{DooerEtAl2004} the authors used a mathematical modeling approach, specif­ically concentrating on rates of change, to replace remedial mathematical approach in summer prep program for remedial mathematics and downstream courses saw an immediate improvement. 

Addi­tionally, we owe our students more than just teaching them algorithmic or even theoretic mathematical material.

In a paper for teachers of American history, Lendol Calder forcefully supports our approach, by quoting the distinguished professor of history, Charles G. Sellers (UC Berkeley):

\begin{quote}
The notion that students must first be given facts and then at some distant time in the future will 'think' about them is both a cover-up and a perversion of pedagogy. \ldots  One does not collect facts he does not need, hang on to them, and then stumble across the
propitious moments to use them. One is first perplexed by a problem and then makes use of the facts to achieve a solution.\cite{Calder2006}
\end{quote}

Put modeling front and center and most importantly put modeling in the driver's seat.
Teacher's are, by the very profession, learners. Indeed, we always say to each other, ``I learn the material best when I teach it.'' As such, learning new mathematics and applications of mathematics present a professional development opportunity. When teaching modeling a teacher has to go along for a ride (sometimes bumpy, but always reward-packed) as students are in the lead all too often and they can go to some new and strange lands, pack some unfamiliar tool boxes, and take some strange transportation modes. It is tiring going in many directions in a class where students are building models. Teachers have to think on their feet, even on their toes! Seeking out modeling material \cite{Winkel2013} offers many rewards, among them the actual finding of useful material for teaching which can enrich students' learning. In that spirit we offer analyses of  several models for an epidemic with full model assessment in light of data.
}

\section*{Epidemic Modeling}

We discuss the traditional SIR model for epidemics and offer a modification, based on Michaelis-Menten kinetics model, to reflect the saturation of infected rates. In both cases we suggest three strategies for estimating parameters on our way to fully assessing the value of each model and each approach. Finally, we compare our results and ask questions.

\subsection*{Weakness in traditional model}

In (\ref{eq:SIR1}) the term $a\cdot S(t)\cdot I(t)$ leads to the following consequences. If we have $I(t) = 1$ infected individuals at time $t$ then the first differential equation for $S(t)$ in (\ref{eq:SIR1}) suggests that the susceptible pool decreases at a rate proportional $S(t)$, i.e. , $-a\cdot S(t)\cdot I(t) = -a\cdot S(t)\cdot 1 = -a\cdot S(t)$, which means that this one individual will infect the same percentage per unit time no matter how many susceptibles there are available, namely $a$, of the susceptibles, $S(t)$. Thus, suppose $a = 0.1$. Then for $S(t) = 1,000$ we would have $(0.1)\cdot 1000 = 100$ new infecteds per unit time while for $S(t) = 100,000$ we would have $(0.1)\cdot 100000 = 1000$ new infecteds per unit time. Indeed, if $S(t)$ were very large then this one individual, $I(t) = 1$, could infect the same percentage of this large population of susceptibles, $S(t)$, {\it no matter how large $S(t)$ is\/}. This is highly unrealistic, for this one infected individual could not contact or transmit  infections to such a large number of susceptibles per unit time.

\subsection*{Strength in modified model with saturation}

In the differential equation for $S'(t)$ in  (\ref{eq:SIR2}) we attempt to be more realistic by replacing the term $a\cdot S(t)\cdot I(t)$ with the term $a\frac{ S(t)}{K_m+ S(t)}I(t)$. The term $K_m+ S(t)$ in the predator-prey literature is an example of a {\it functional response\/} \cite{Holling1959} or rate of prey intake by the predator as a function of the prey levels. This functional response suggests that as $S(t)$ gets larger and larger there is a limiting value of the term
representing the rate of infection, indeed, $\lim_{t->\infty} \frac{ S(t)}{K_m+ S(t)}  = 1$. College Joe (predator) can consume
only so many pizzas (prey) per unit time even if presented the whole pizza store's supply! In our situation we do not have predator or prey, but rather infected and susceptible, respectively.

If we consider the full term $a\frac{ S(t)}{K_m+ S(t)}I(t)$  in the first differential equation for $S'(t)$  in (\ref{eq:SIR2}) and let
$S(t)$ get larger and larger, i.e.~the available susceptible pool level increases indefinitely, we see that
$a\frac{ S(t)}{K_m+ S(t)}I(t) \to a \cdot I(t)$, which means that the rate of infection or decrease in susceptibles, depends only
on the number of infected individuals. Indeed, $a$ can be considered as the maximum infection rate of infecting susceptible individuals attainable per individual infected member of the population. The number $K_m$ signifies the number of susceptibles available to the infected population at which infection takes place at exactly half this maximum rate, $a$. We illustrate this in Figure~\ref{fig:6-019-EpidemicArticleHolling} in which $K_m = 40$. We see that when the susceptible population level is at $S = 40$ then the functional response level or maximum infection rate is at 50\% or 0.50.

\section*{Parameter Estimation Efforts}

Now that we have two models for the spread of the disease, (\ref{eq:SIR1})  and (\ref{eq:SIR2}) , we seek to ascertain if either or both render good models for the boarding school data in Table~\ref{table:data} by estimating the parameters for both and fitting the resulting models to the data in Table~\ref{table:data}. We proceed in three ways,  based on minimization of the sum of square errors between the model prediction of infected individuals and the observed number of infected individuals over the entire 14 day period of the epidemic.

First, in the case of (\ref{eq:SIR1}) we consider the parameter surface, ${\rm ParSpace} = (a, b, SSE(a b))$ in three dimensional space. formed by using the parameters $a$ and $b$ to numerically solve the nonlinear system of differential equations in (\ref{eq:SIR1}). We then compute the sum of square errors, $SSE(a, b)$ between the model prediction of infected individuals and the observed number of infected individuals over the entire 14 day period of the epidemic for the parameter values of $a$ and $b$. The totality of these points, $(a, b, SSE(a, b))$, form our sum of square errors function or parameter surface in this parameter space. We seek to determine the parameters $a$ and $b$ which minimize $SSE(a, b)$.

We can do this minimization in several ways.

\begin{enumerate}
\item 
Select reasonable initial estimates of $a$ and $b$ from a range of values, say $a_{\rm min} \leq  	a \leq  a_{\rm max}$ 
and $b_{\rm min} \leq  	b \leq  b_{\rm max}$. Then use trial and error by numerically obtaining solutions of (\ref{eq:SIR1}) and comparing the resulting model plot with the data plot as seen in Figure~\ref{fig:6-019-EpidemicArticleEarlyEstimate}.


\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=11.5cm]{6-019-EpidemicArticleEarlyEstimate}



\parbox{12. true cm}{\caption{\label{fig:6-019-EpidemicArticleEarlyEstimate} An early estimate for the model in (\ref{eq:SIR1}) using $a = 0.005$ and $b = 0.5$, resulting in a sum of quares error of $SSE(0.005, 0.5) = 318,089$. This is not a particularly good estimate!}}
\vspace*{-.3 true cm}
\end{figure}

\item We can compute points on the three dimensional surface, ${\cal S}$ defined by (\ref{eq:CalS}):


\begin{equation}\label{eq:CalS}
{\cal S}  = \{(a, b, SSE(a, b)\vert a_{\rm min} \leq  	a \leq  a_{\rm max}, \,\, b_{\rm min} \leq  	b \leq  b_{\rm max}\}\, .
\end{equation}

We do this over a grid with a small mesh size, perhaps with thousands of points and hence thousands of solutions of our differential equations system (\ref{eq:SIR1}) for each paramter pair $(a, b)$, along with an accompanying computation of the sum of square errors betweeen model values and data values. From this set of values, ${\cal S}$, we extract the triple $(a, b, SSE(a, b))$ whose $SSE(a, b)$ is smallest and offer up the corresponding parameters $a$ and $b$ as our ``best'' model parameters. We can then compare our model using these optimal parameter values to our data for a visual confirmation of our modeling process.
\item We can use a gradient approach of steepest descent on the surface, ${\cal S}$, to determine values of $a$ and $b$ which render a minimum $SSE(a, b)$. We compute the gradient by estimating the respective partial derivatives in each direction, $a$ and $b$, using symmetric differences, and forming a unit vector in the opposite direction of the resulting unit gradient vector, for the gradient vector
points in the direction of steepest ascent and we want to descend to a minimum; hence we go in the opposite direction, i.e.~{\it steepest descent}
 
We  compute the partial derivatives numerically using a symmetric difference approach, in (\ref{eq:part}). For example, we compute $\frac{\partial SSE(a,b)}{\partial a}$  as:

\begin{equation}\label{eq:part} 
\frac{\frac{SSE(a+\Delta a,b)-SSE(a,b)}{\Delta a} + \frac{SSE(a,b)-SSE(a-\Delta a,b)}{\Delta a}}{2} 	= \frac{SSE(a + \Delta a, b) - SSE(a - \Delta a, b)}{2 \Delta a}  \, .
\end{equation}


It should be pointed out that each evaluation of $SSE(a, b)$ involves placing the parameters $a$ and $b$ into the system of differential equations (\ref{eq:SIR1}), solving the system numerically, and computing the sum of square errors between the resulting model for the number of infected individuals with this particular parameter set $a$ and $b$ and the actual data values for the number of infected individuals. Thus, each step is computational intense. Moreover,we may change the step size as we near the minimum value of $SSE(a, b)$.

We can also do something comparable for a grid and a gradient search for the saturation model, (\ref{eq:SIR2}), except that we are then in a four-dimensional world and our ``surface'' is defined by (\ref{eq:CalSK}):

\begin{equation}\label{eq:CalSK}
{\cal SK}  = \{(a, b, SSE(a, b)\vert a_{\rm min} \leq  	a \leq  a_{\rm max}, \,\, b_{\rm min} \leq  	b \leq  b_{\rm max}, \,\, K_{m, {\rm min}} \leq  K_m \leq  K_{m,{\rm max}}\}\, .
\end{equation}

Nevertheless, (1) we can adopt a trial and error approach - more computational choices now as there are three parameters $a$, $b$, and $K_m$, but we may have a sense of their value and meaning to proceed intelligently, not just blindly; (2) proceed to form a three dimensional grid as in the first model with only two paramteres, $a$ and $b$,  bu in this model three parameters, $a$, $b$, and $K_m$, and search for a minimum $SSE(a, b, K_m)$ in this case; and (3) compute gradients and use the method of steepest descent to move down the surface, ${\cal SK}$ (in four dimensional space!), efficiently to a minimum $SSE$ value and hence a best fitting set of parameters for our model in (2).
\end{enumerate}



\section*{ACTIVITY}

For both models (\ref{eq:SIR1}) and (\ref{eq:SIR2}) determine the best model using one or more of the strategies:
\begin{enumerate}
\item Trial and error approach in an attempt to guess optimal parameter values for modeling the data.
\item  Form a  grid  with  a reasonable range for the parameters   and search for a minimum $SSE$ as a function of these parameter values.
\item Compute gradients and use the method of steepest descent to move down the surface, whose coordinates are the parameters and the sum of squared errrors, efficiently to a minimum $SSE$ value.
\end{enumerate}



\teacher{
\section*{NOTES FOR TEACHERS}

As supplemental materials for this article we offer two Mathematica notebooks, 6-019-T-Mma-EnablingEpidemicExploration-MM-TeacherVersion and  6-019-T-Mma-EnablingEpidemicExplora- \linebreak tion-NonMM-TeacherVersion each with complete analyses. We also offer pdf versions of each file.

Here we give results on the various methods (1) trial and errors, (2) grid, and (3) gradient, with discussion of resolution and stopping rules for both (a) SIR model and (b) SIR model with Michaelis-Meriten saturation terms. We present graphics in Figure~\ref{fig:6-019-EnablingEpidemicContourFigure}  to show convergence in two dimensions over contour plot.

\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=11.5cm]{6-019-EnablingEpidemicContourFigure}



\parbox{12.5 true cm}{\caption{\label{fig:6-019-EnablingEpidemicContourFigure} 
Contour plot of SIR model without Michaelis-Meuten terms (1) in which the results $(a, b, SSE(a, b))$ of the grid search are depicted as contour levels and the gradient search path is laid over the contour plot. The point of convergence, $(a = 0.00258$,  $b = 0.484$, $SSE(a, b,) = 5129.79$) (see Table~\ref{table:summary}) is realized by both methods.}}

\end{figure}

The mesh for the grid approach is over intervals for each parameter that are reasonably close to estimate final values due to experimentation and in each direction we take 200 intervals. Thus in the case of $a$ and $b$ in the SIR model we have some 40,000 grid values of $(a, b, SSE(a, b))$ while in the Michaelis-Menten approach we take 20 intervals in each of the three parameters $a$, $b$, and $K_m$, for 8,000 values of $(a, b, K_m, SSE(a, b, K_m))$.

For the gradient approach we usually take about 1,000 iterations with the step size, being cut by a factor of 0.99 at each iteration. Thus the step size gets so small that there is hardly any difference in our $SSE$ values after a few hundred iterations. This is more than sufficient to get us very close to our optimal values. Quite obviously, since the grid approach is a function of the number of intervals and the original interval width we would expect this inflexibility to not give us as good a value as the gradient which is more dynamic, ever seeking to better itself.

Figure~\ref{fig:6-019-EnablingEpidemicContourFigure} shows the result of both the grid and gradient method which converge to the same values of $a$ and $b$ in the SIR model. This is reassuring.

\begin{table}
\begin{center}
\begin{tabular}{|l|c|c|c|c|} \hline \hline
\multicolumn{1}{|c|}{Model/Method}	 &&&&\\ \hline\hline	   
SIR Model	&a	&b	&-	&$SSE(a, b)$	\\ \hline    
Grid	&0.00258	&0.484&	-	&5181.43	\\ \hline  	   
Gradient	&0.00257	&0.483&	-	&5129.79	\\ \hline   \hline 	   
SIR Model with Michaelis-Mentert	&$a$	&$b$	& $K_m$ & $SSE(a, b, K_m)$		\\ \hline     
Grid&	5.875	&0.51	&1600.00	&10245.8		\\ \hline     
Gradient&	8.09983	&0.50046	&2467.02&	8313.66	 	\\ \hline   \hline
\end{tabular}


\parbox{11.8 true cm}{\caption{\label{table:summary}Summary of both Model (SIR with and without Michaelis-Menten terms) and Method (grid and gradient) with parameters and sum of square error terms.}}
\end{center}

\end{table}




\section*{Issues}

What does Table~\ref{table:summary} offer us? When we examine the resulting best model using the gradient approach in both cases, without and with the Michaelis-Menten terms in Figure~\ref{fig:6-019-EpidemicModelDataPlotNonMichaelisMenten} and Figure~\ref{fig:6-019-EpidemicModelDataPlotWithMichaelisMenten}, respectively, we see that without Michaelis-Menten terms our plot looks like it holds to the data better. However, that really is an aesthetic call and upon numerical comparisons we see that our $SSE$ in the case of no Michaelis-Menten terms ($SSE(a, b,) = 5129.79$) is better than that in the case of using Michaelis-Menten terms ($SSE(a, b, K_m) = 8313.66$). Now which is the better model? Well, perhaps pictures and numbers cannot really tell us here.

We just might need to decide on which is the better model on our understanding of the models themselves. In both models SIR model without Michaelis-Menten terms (1) and the SIR model with Michaelis-Menten terms (2) the coefficient $b$ has the same significance and units, namely, $b$ is the exponential decay rate for the Infected individuals  moving into the Recovered individual compartment. Both models give about the same value for $b$ (which is reassuring), namely $b = 0.484$ for (1) and $b = 0.50046$ for (2). Thus it is in an examination of the terms representing the effect of Infected Individuals on Susceptible individuals that we must concentrate.

\subsection*{Significance of Parameters}

The units on the respective parameters are important to consider and we review them in Table~\ref{table:best}.

\section*{Comparing Models}

We have two models for our epidemic data and both appear reasoned and have close sum of square errors. The first, (\ref{eq:SIR1}), has two parameters $a$ and $b$ while the second, (\ref{eq:SIR2}), has three parameters, $a$, $b$, and $K_m$.

The first has a sum of square errors of 5129.79 while the second has a sum of square errors of
8313.66. Some texts refer to the sum of square errors ($SEE$) as the residual sum of errors ($RSE$). Now certainly if we increase the number of parameters we can get a better fitting model in terms of reducing the $SSE$. Indeed, if we took a $14^{\rm th}$  degree polynomial we could fit it to go through all the 14 data points exactly and thus our $SSE$ would be 0. However, this model would not be valuable in any predictive sense, nor would it add to our understanding of the phenomenon we are modeling. Moreover, it would not explain anything as to the actions of the underlying model. So there is a trade-off, namely, if we increase the number of parameters in our model we can decrease our $SSE$, but to what avail in our understanding of the meaning of our model. 






We can ask which model is best, in general. Namely, which model keeps both the $SSE$ and the number of parameters low.
In 1974 a formula was published \cite{Akaike1974} which would allow us to use both $SSE$ and number of parameters at the same time to compare models. The statistic is called the {\it Akaike Information Criterion\/} or $AIC$. A reasonable formulation (\ref{eq:AIC}) for our purpose is given in \cite[p.~138]{Akaike1974}:

\begin{equation}\label{eq:AIC}
AIC(n,k,SSE) = n\cdot \ln\left(\frac{SSE}{n}\right) + 2(k + 1)\, ,
\end{equation}

\noindent where $n$ is the number of observations, $k$ is the number of parameters to be estimated, and $SSE$ is the least sum of square errors.
$AIC$ is used to compare models and the model with the smallest $AIC$ is considered best of the set. $AIC$ does not tell us if the model is good, only if it is the best of the ones being considered.

\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=11.5cm]{6-019-EpidemicModelDataPlotNonMichaelisMenten}



\parbox{12.5 true cm}{\caption{\label{fig:6-019-EpidemicModelDataPlotNonMichaelisMenten} Plot of SIR model without Michaelis-Menten terms (\ref{eq:SIR1}) in which the resulting parameters $a = 0.00258$ and $b = 0.484$ give a minimum $SSE(a, b,) = 5129.79$.}}

\end{figure}

 
 \begin{table}[!h]
 \begin{center}
 \begin{tabular}{|c|c|c|c|c|}\hline
Type &	$a$&	$b$&	$K_m$&	$SSE$	   \\ \hline
Non-Michaelis-Menten	&$\frac{1}{{\rm 	time}\cdot {\rm Infecteds}}$		&$\frac{{\rm Recovereds}}{{\rm 	time}\cdot {\rm Infecteds}}$	& - &  \\ 
 & 0.00258 & 0.484 && 5,129.79 \\ \hline	 	   
Michaelis-Menten	&$\frac{{\rm Susceptibles}}{{\rm 	time}\cdot {\rm Infecteds}}$		&$\frac{{\rm Recovereds}}{{\rm 	time}\cdot {\rm Infecteds}}$	& Susceptibles  & \\ 
&8.09983 & 0.50046& 2,4667.02 & 8,313.66 \\ \hline
\end{tabular}

\caption{\label{table:best} Units and best values of parameters in the two models under consideration.}
\end{center}
\end{table}




We offer Table~\ref{table:AICStuff} for the $AIC$ values for our two models, the non-Michaelis-Menten and the Michaelis-Menten model. From Table~\ref{table:AICStuff} we see that even when we introduce more parameters and saturation for infecting term into the standard SIR model it does not improve the original model.

 
 \begin{table}[!h]
\begin{center}
\begin{tabular}{|c|c|c|c|c|} \hline
Type	&$n$&		$k$	&$SSE$	&$AIC(n,k,SSE)$	\\ \hline   
Non-Michaelis-Menten	&14	&2	&5,518.13	&89.6743	\\ \hline     
Michaelis-Menten&	14&	3	&8,313.66&	97.4124	 \\ \hline  
\end{tabular}

\caption{\label{table:AICStuff} Summary of $AIC$ values for SIR without and with Michaelis-Menten terms.}
\end{center}
\end{table}

\begin{figure}[h]
\centering

\vspace*{-.1 true cm}
\includegraphics[width=11.5cm]{6-019-EpidemicModelDataPlotWithMichaelisMenten}



\parbox{12. true cm}{\caption{\label{fig:6-019-EpidemicModelDataPlotWithMichaelisMenten} Plot of SIR model with Michaelis-Menten terms (\ref{eq:SIR2}) in which the resulting parameters $a = 8.09983$, $b = 0.50046$, and $K_m = 2467.02$, give a minimum $SSE(a, b, K_m) = 8313.66$.}}

\end{figure}


Our conclusion from these analyses is that while the introduction of the Michaelis-Menten saturation term does not statistically permit us to say the newer model in this case is better we have every reason to accept the plausibility of the newer model because of its saturation nature and how it makes more sense than the non-Michaelis Menten model.


}





\begin{thebibliography}{99}

\teacher{
\bibitem{Akaike1974} Akaike, H. 1974. A new look at the statistical model identification. {\it IEEE Transactions on Auto­matic Control\/}. 19(6): 716-723.
}

\teacher{
\bibitem{Calder2006} Calder, L. 2006. Uncoverage: Toward a Signature Pedagogy for the History Survey. {\it The Journal of American History\/}. 92(4): 1358-1370.
}

\bibitem{CDC1978}
Communicable Disease Surveillance Center. 1978. News and Notes: Influenza in a Board­ing School. {\it British Medical Journal\/}. 1(6112): 586-590. \href{https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1603269/pdf/brmedj00115-0064.pdf}{https://www.ncbi.nlm.nih.gov/pmc/articles/PMC1603269/pdf/brmedj00115-0064.pdf}. 

\teacher{
\bibitem{DooerEtAl2004} Dooer, H~.M., J.~B.~Arleback, and A.~C.~Staniec. 2014. Design and Effectiveness of Modeling-Based Mathematics in a Summer Bridge Program. {\it Journal of Engineering Education\/}. 103(1): 92-114.
}



\bibitem{Holling1959} Holling, C. S. 1959. The components of predation as revealed by a study of small-mammal pre­dation of the European pine sawfly. {\it The Canadian Entomologist\/}. 91(5): 293-320.

\teacher{
\bibitem{Ledder2013} Ledder, G. 2013. {\it Mathematics for the Life Sciences\/}. New York: Springer.
}

\bibitem{MillerHelms2010} Miller, S. and J.~Helms. 2010. The Disease that Infected Half the World: Mathematics, Biology, and History. {\it PRIMUS - Problems, Resources, and Issues in Mathematics Undergraduate Studies\/}. 20(3): 245-260. 


\bibitem{Miller2014} Miller, S. 2014. Modeling a Nonlethal Influenza Epidemic. 6-001-S-Epidemic. 

\bibitem{Murray1993} Murray, J. 1993. {\it Mathematical Biology, Second Corrected Edition\/}. New York: Springer-Verlag.

\bibitem{Switkes2003} Switkes, J. 2003. A modified discrete SIR model. {\it The College Mathematics Journal\/}. 34(5): 399-­402.

\teacher{
\bibitem{Winkel2011} Winkel, B. 2011. Parameter Estimates in Differential Equation Models for Population Growth. {\it PRIMUS\/}. 21(2): 101-129.


\bibitem{Winkel2013} Winkel, B. 2013. Browsing Your Way to Better Teaching. {\it PRIMUS\/}.
 23(3): 274-290.
}

\teacher{
\bibitem{Winkel2015} Winkel, B. J. 2015. 6-001-T-EXCEL-Epidemic-TeacherVersion. 

\bibitem{Winkel2018} Winkel, B. 2018. Parameter Estimation Analysis for Epidemic Models. Supplemental Mathematica Modeling Scenario for this article found at the SIMIODE website.
}
\end{thebibliography}






\end{document}










v