A Model for the Spread of Infectious Diseases with Application to COVID-19

: Given the present pandemic caused by the severe acute respiratory syndrome coronavirus 2 or SARS-CoV-2 virus, the authors tried ﬁtting existing models for the daily loss of lives. Based on data reported by Worldometers on the initial stages (ﬁrst wave) of the pandemic for countries acquiring the disease, the authors observed that the logarithmic rendering of their data hinted the response of a ﬁrst-order process to a step function input, which may be modeled by a three-parameters function, as described in this paper. This model was compared against other similar, log(N)-class of models that are non-compartmental type (such as the susceptible, infected, and removed, or SIR models), obtaining good ﬁt and statistical comparison results, where N denotes the cumulative number of daily presumed deaths. This simple ﬁrst-order response model can also be applied to bacterial and other biological growth phenomena. Here we describe the model, the numerical methods utilized for its application to actual pandemic data, and the statistical comparisons with other models which shows that our simple model is comparatively outstanding, given its simplicity. While researching the models available, the authors found other functions that can also be applied, with extra parameters, to be described in follow-on articles.


Introduction
Population growth has been studied in many fields, mostly in biology. Initially studied by Robert Malthus on population growth, this initial model has been transformed and adapted to particular cases which are commonly referred to as logistic equations, described by sigmoidal-type functions that commence with a small growth rate but which increases until reaching a maximum rate, and after this peak in growth rate, it diminishes and the curve flattens out asymptotically to a rate near zero, where the population reaches some maximum value, called carrying capacity, horizontal asymptote, saturation value, etc. The models that describe these growth phenomena are described in numerous papers and classic books on the subject like [1,2] and others.
One of the main reasons the authors chose of using a simple growth model (which is generally applied to bacterial growth) arises from the uncertainty in the reported numbers of infected, recovered, and removed (recovered and/or dead) people during the onset of the COVID-19 pandemic. There are many recently published papers about mathematical modeling and analysis of COVID-19 [3][4][5][6][7][8][9][10]. The model presented here was derived based on data collected at the onset (or "initial-wave") of the pandemic, and our study, restricted to this initial-wave, does not include the subsequent waves (as are being experienced presently), which would require a highly complex dynamical model to track.
Several of the countries for which we tried fitting our model have reset their numbers, sometimes on a weekly basis, and in other cases, a long time-interval of data collected The relaxation of precautionary measures in many countries has produced a deviation of the expected or "well-behaved" rate of presumed deaths curve after its peak has been reached. The tail of the decaying rate curve, instead of following a monotonic decay, has changed (in many countries) its curvature again. For example, in the USA, the peak of death rates occurred around the third week of April 2020, and started to decay, until a change in curvature (using the IHME model prediction data) occurred around 28 June 2020. Because of these changes being experienced in many countries, the authors decided to utilize only data from the pandemic onset (data collected since 21 January 2020) up to 31 May 2020. Due to the variations in the expected (or well-behaved) response, the modeling of such variations would require a dynamical update model such as the one utilized by the IHME, with an inherent increase in model complexity. Our model in contrast, however, fitted the onset of the pandemic quite accurately, as compared to other similar log(N) models and the IHME model in the early onset of the pandemic.
The outline of the paper is as follows. Our mathematical model based on a first-order process to a step function input is put forward in Section 2. In Section 3 our model is compared to data reported by Worldometers [11] and compared against other similar, log(N)-class of models. In Section 4, we point out the importance of different model approaches. Concluding remarks form Section 5.

Mathematical Model
The starting point of our model is the observation that the logarithm of the actual recorded data by [11] on the number of deaths follows a curve that can be associated with a first order process, that obtained from a first order differential equation, as the response of a first order system to a Heaviside (or step) function type of excitation. Figure 1 describes, pictorially, the steps taken to construct the model. Figure 1a depicts the shape of the curve (sigmoid) of the actual reported data, y a (t), which is the cumulative number of deceased people as a function of elapsed time, in days, since the outbreak of the pandemic for a specific country. Figure 1b depicts the natural logarithm of the actual recorded data obtained by a process P 1 , such as y 1 (t) = ln(y a (t)), and Figure 1c depicts the linearization of the curve in Figure 1b, i.e., y 2 (t) =â 0 +â 1 t, obtained by a linearization process P 2 . tion of the expected or "well-behaved" rate of presumed deaths curve after its peak has been reached. The tail of the decaying rate curve, instead of following a monotonic decay, has changed (in many countries) its curvature again. For example, in the USA, the peak of death rates occurred around the third week of April 2020, and started to decay, until a change in curvature (using the IHME model prediction data) occurred around 28 June 2020. Because of these changes being experienced in many countries, the authors decided to utilize only data from the pandemic onset (data collected since 21 January 2020) up to 31 May 2020. Due to the variations in the expected (or well-behaved) response, the modeling of such variations would require a dynamical update model such as the one utilized by the IHME, with an inherent increase in model complexity. Our model in contrast, however, fitted the onset of the pandemic quite accurately, as compared to other similar log(N) models and the IHME model in the early onset of the pandemic.
The outline of the paper is as follows. Our mathematical model based on a first-order process to a step function input is put forward in Section 2. In Section 3 our model is compared to data reported by Worldometers [11] and compared against other similar, log(N)-class of models. In Section 4, we point out the importance of different model approaches. Concluding remarks form Section 5.

Mathematical Model
The starting point of our model is the observation that the logarithm of the actual recorded data by [11] on the number of deaths follows a curve that can be associated with a first order process, that obtained from a first order differential equation, as the response of a first order system to a Heaviside (or step) function type of excitation. Figure 1 describes, pictorially, the steps taken to construct the model. Figure 1a depicts the shape of the curve (sigmoid) of the actual reported data, () a y t, which is the cumulative number of deceased people as a function of elapsed time, in days, since the outbreak of the pandemic for a specific country. Figure 1b   The curves shown in Figure 1 are mathematical idealizations. These curves are actually described by the discrete daily values, which exhibit daily variations around these mathematically idealized curves or functions. To obtain the idealized curves, the data in the linearized curve in Figure 1c is approximated, using a weighted least squares fitting process, using a line equation, and then, the inverse processes  The curves shown in Figure 1 are mathematical idealizations. These curves are actually described by the discrete daily values, which exhibit daily variations around these mathematically idealized curves or functions. To obtain the idealized curves, the data in the linearized curve in Figure 1c is approximated, using a weighted least squares fitting process, using a line equation, and then, the inverse processes P −1 2 and P −1 1 are applied successively to obtain the idealized first-order model (curve in Figure 1b), and the "actual" function ( Figure 1a) that fits the actual recorded data. Some key important information about the pandemic can be obtained from applying the processes shown in Figure 1, and by their corresponding inverse fitting processes: (1) An approximating function based on the proposed for the actual data that enables to predict near future outcomes and also enables to backtrack to estimate the initial values that might have been in those cases where the early or initial data on the pandemic outbreak is not available or was not reported or recorded. (2) A time-constant ( τ) of the first order system that determines the time for which the number of deceased people will reach a particular percentage of the peak value. Recall that the response of a first-order system would reach a plateau value asymptotically. For example, the population will reach 99.33% of the plateau (or saturation) value y ∞ = y 1 (∞) for a period of time equal to five time constants or t = 5 τ; 99.75% of y ∞ for t = 6 τ, etc. (3) The asymptote's (maximum) value that corresponds to the maximum number of people that would be deceased in the long run or for a long period of time. (4) The inflection point of the estimated curveŷ a (t) at which the sigmoid changes concavity in Figure 1a. This point corresponds also to the time for which the maximum rate of deaths per day (i.e., dŷ a (t)/d t) occur. We next describe the two processes P 1 and P 2 that convert the actual reported data of Figure 1a into that in Figure 1b (our model), and that in Figure 1b into the linear form in Figure 1c, respectively.
Process P 1 : This is a straightforward process that consists simply in taking the natural logarithm of the actual reported data y a (t) and creating a new array with the natural logarithm of the actual data values versus time. That is, one generates the array of values Since this data resembles the response of a first order-system to a Heaviside (or step) input function u(t) or stepwise excitation given by y ∞ u(t), where y ∞ is the horizontal asymptote or saturation value, the curve in Figure 1c is given by where the unit-step or Heaviside function u(t) is defined piecewise as u(t) = 1, t ≥ 0 0, t < 0 .
In (2), t represents time in units of [day], α is a growth rate in units of [1/day], and its reciprocal τ = α −1 (days), the time-constant corresponding to the first order-process or system with initial condition y 1 (0) = 0. This implies that y 1 (0) = ln(y a (0)) = 0 which requires that y a (0) = 1. This indicates that our coordinates require to start at the point where the first death occurs, that is, our first point is (t, y a ) = (0, 1). Therefore, the time t elapsed is measured from the day the first death occurred.
Process P 2 : This process linearizes the first order response curve, by transforming Equation (2) by means of a simple algebraic manipulation and a change of variables as follows: For t ≥ 0, one has Letting one can obtain the linear equationŷ 2 (t) = αt, after determining α via least squares process, which as shown in Figure 1c. A more general model (3-parameters: y ∞ , a, t d , where y ∞ is the asymptote value, a is the shape parameter, and t d is a time delay) for a first order process response would include an initial value y i = y 1 (t d ), where t d is a time-delay parameter of the process relative to some zero time reference and where the first data point exists, at (t, y 1 ) = (t d , y i ). The shifted version of the first order model shown in Figure 2b can be written in a more general form as FOR PEER REVIEW 6 of 20 From (7), after the line coefficients have been determined via weighted least squares, one can write Exponentiating both sides of the second equality in (8) and solving for In (9) We note that in this study we use weighted least squares, where the weights get larger towards the later data and smaller towards the initial data values. The hypothesis here is that the residuals of the natural logarithm of the data and the fitting function become smaller as the rate of increase becomes smaller, towards the horizontal asymptote, so we weight heavier towards the region of smaller residuals. We tried several simple weighing functions that minimize the sum of squares of the residuals between the actual and the fitted data, such a linear, quadratic, and higher degree weighting functions, where the sum of the weights in each case is normalized to one. That is, n Figure 2. A more general model: (a) Actual recorded data y a (t) and its approximating functionŷ a (t). (b) Logarithm of the sigmoidal curve and its approximating functionŷ 1 (t). (c) Linearization of the y 1 (t) curve and its approximating functionŷ 2 (t).
This function is essentially the response of a first-order system to a Heaviside [or step u(t)] input function. Derivation of this response is explained in the Supplementary Materials section, under Mathematics of first-order differential equation. Since this system is linear-time-invariant (LTI), any input delayed in time would correspond to the same delayed response. That is, one gets the response of a non-anticipative, causal system. More detailed description of LTI systems is given in reference [20].
The linearization (process P 2 ) for this more general model is obtained algebraically as follows: From (5), for t ≥ t d , one has ln In (6), the time t * = t − t d is a time shifted time. Letting one can obtain a linear least-squares approximation equationŷ 2 (t * ) = αt * = α(t − t d ) = a 1 t +â 0 shown in Figure 2c, The right pointing arrows in Figure 2 indicate the two processes that transform the actual data into the curve corresponding to our model and the linearization of the data, while the left pointing arrows indicate the inverse processes, after the linearized data has been fitted via a weighted least Squares process.
After the coefficient estimateâ 0 andâ 1 of the line equationŷ 2 (t) have been determined via a weighted least squares process, then the inverse processes P −1 2 and P −1 1 are applied to obtain the approximating or fitting functionsŷ 1 (t) andŷ a (t), respectively. A description of the inverse processes follows.
From (7), after the line coefficients have been determined via weighted least squares, one can writeŷ Exponentiating both sides of the second equality in (8) and solving forŷ 1 (t) yieldŝ Challenges 2021, 12, 3 6 of 19 In (9) the valuesâ 1 = α,â 0 = −α t d have been replaced back, andŷ i is the initial valuê y 1 (t = t d ). Once the functionŷ 1 (t) has been determined, thenŷ a (t) is obtained, simply, aŝ y a (t) = exp(ŷ 1 (t)) (10) We note that in this study we use weighted least squares, where the weights get larger towards the later data and smaller towards the initial data values. The hypothesis here is that the residuals of the natural logarithm of the data and the fitting function become smaller as the rate of increase becomes smaller, towards the horizontal asymptote, so we weight heavier towards the region of smaller residuals. We tried several simple weighing functions that minimize the sum of squares of the residuals between the actual and the fitted data, such a linear, quadratic, and higher degree weighting functions, where the sum of the weights in each case is normalized to one. That is, Linear weighting function: For meeting requirement given by Equation (11), the linear weights, as a function of the data point, are given by Quadratic weighting function: For meeting requirement given by Equation (11), the quadratic weights, as a function of the data point, are given by Cubic weighting function: For meeting requirement given by Equation (11), the cubic weights, as a function of the data point, are given by Since all weights are positive, the weighted linear least squares parameter estimateŝ a 0 andâ 1 are obtained by solving the normal equations, that is, In (15), W is the diagonal matrix: In (15), a = â 0 a 1 is the vector of coefficients, . . .
is the vector of transformed (linearized) values y 2 (t) obtained via (4). We note that the value for y ∞ is obtained by a least squares approach as well, by minimizing the square of the residuals r i (t) = y i (t) −ŷ i (t), that is, by finding a value y ∞ that minimizes RSS = r T r where r is a column vector of the residuals and T denotes transposition. Additional details on the fitting process such as error estimates for the parameters and the corresponding sensitivities are provided in Appendix A.

Death Rates
The death rate in number of deaths per day is calculated by taking the derivative of the estimated actual accumulated death data computed by (10). Taking the time derivative of (10) yields In (17), the time derivative ofŷ 1 (t) is obtained from (9) as The graph of (17) shows an initial gradual increase in the values of the rate or time derivative in units of [deaths/day], then it reaches a peak or maximum value corresponding to the inflection point of the sigmoid and after this point commences to decrease towards the value of zero, asymptotically. The graphs of the death rates are also given for each of the selected countries in the Supplementary Materials section and for Germany in Figure 3 in the Results section. These rate graphs are useful to find the time when the peak rate occurs, date from which the effects of pandemic begin to diminish. The graph of (17) shows an initial gradual increase in the values of the rate or time derivative in units of [deaths/day], then it reaches a peak or maximum value corresponding to the inflection point of the sigmoid and after this point commences to decrease towards the value of zero, asymptotically. The graphs of the death rates are also given for each of the selected countries in the Supplementary Materials section and for Germany in Figure 3 in the Results section. These rate graphs are useful to find the time when the peak rate occurs, date from which the effects of pandemic begin to diminish.

Results
The model was applied to the actual data that was collected on a daily basis from the web site given in reference [11] which provides historical and daily updated values of the COVID-19 infections and presumed deaths for numerous countries. Initially, we collected the data for a limited set of countries, those that were reported initially, starting with China. The number of countries that were included in our data set increased in time as the pandemic spread over the world. In Section 3.1, our proposed First Order model is fitted to each of the countries on the initial data set. In Section 3.2, already published models are fitted to this data set, and the results are compared to those in our proposed first order model.

Results
The model was applied to the actual data that was collected on a daily basis from the web site given in reference [11] which provides historical and daily updated values of the COVID-19 infections and presumed deaths for numerous countries. Initially, we collected the data for a limited set of countries, those that were reported initially, starting with China. The number of countries that were included in our data set increased in time as the pandemic spread over the world. In Section 3.1, our proposed First Order model is fitted to each of the countries on the initial data set. In Section 3.2, already published models are fitted to this data set, and the results are compared to those in our proposed first order model.

Simulation Results
The data we collected for this article span from 21 January 2020 (earliest data for China) to 31 May 2020, which we considered sufficient, in number of data points, to test our model within varying intervals of time for verifying the model estimates and the characterization of the model parameters with the values of the reported data [11].
Daily data records for ten of the countries utilized in the simulations are given in Table S1 in the Supplementary Materials section. The graphs corresponding to the actual data and the estimated model data [y a (t) andŷ a (t); y 1 (t) andŷ 1 (t); y 2 (t) andŷ 2 (t)]; the horizontal asymptote value y ∞ ; and the Process-P 1 time constant τ = 1/α [days] for the curveŷ 1 (t) for each country listed in Table S1 are also illustrated in the Supplementary Materials section. A fraction (start and end) of Table S1 is repeated here just for reference in Table 1: The graphs that were generated for several of the countries listed in Table S1 of the Supplementary Materials section include the following graphs: (1) Actual data y a (t) and the actual data model estimates curveŷ a (t).
(2) The natural logarithm of the actual data y 1 (t) and the model estimate curveŷ 1 (t).
(3) The linearization of the actual data logarithm y 2 (t) and the linear fitŷ 2 (t). (4) The derivative of the actual data model dŷ a (t)/d t.
(5) Additional graphs of interest, like fit residuals n some of the cases analyzed.
The numerical computations were performed using MATLAB (from MathWorks) scripts written by the authors, which included some high-level (canned-in) commands. A Summary of modeling results/parameters, for the given data recorded in Table 1 for the ten countries selected are provided in Table 2. Very good fits of the proposed First Order model with a RSS value below one were obtained for CN, SK, and the US; good fits with an RSS value between one and two were obtained for IT, IR, SP, DE, UK, and NE. The worst fit was obtained with an RSS value larger than 2 for FR.
For completeness, we provide below a sample of the modeling process for the data corresponding to Germany (DE), corresponding to

Comparison to Other Models
Other log(N) models besides the model presented here, which we call First-Order model given byŷ whereŷ are presented and compared next. The model comparison graphs and the modified forms of these models are also described in Appendix B. These models have been widely applied (and compared) for a diversity of population growth phenomena. The following models are reported in [12,17], and [21] and include: Gompertz model which is a three-parameters curve Richards model which is a four-parameters curve Logistic model, which is another three-parameter curve Stannard model, which is a four-parameter curve and Schnute model another four-parameters curve Parameters like a, b, c, ν, τ, k, l, p, y 1 , τ 1 , and τ 2 in these models denote shape, shift, and/or scale parameters (see first column in Table 3).   Challenges 2021, 12, x FOR PEER REVIEW Figure 6. Data-Rate corresponding to the derivative of the estimated data in Figure 3. The p value occurred on 13 April 2020, according to the model.

Comparison to Other Models
Other log(N) models besides the model presented here, which we call First model given by are presented and compared next. The model comparison graphs and the modified of these models are also described in Appendix B. These models have been widely a (and compared) for a diversity of population growth phenomena. The following m are reported in [12,17], and [21] and include: Gompertz model which is a three-parameters curve Richards model which is a four-parameters curve Logistic model, which is another three-parameter curve Figure 6. Data-Rate corresponding to the derivative of the estimated data in Figure 3. The peak value occurred on 13 April 2020, according to the model. These five models, (21) to (25) and the proposed First-Order model, (19) and (20) were next compared to each other using the logarithm of the original data, y 1 (t) = ln(y a (t)), for the Netherlands. Figure 7 depicts the corresponding First-Order model fit (solid red line) to the data (blue dots).
These five models, (21) to (25) and the proposed First-Order model, (19) and (20) next compared to each other using the logarithm of the original data,  1 ( ) a y t ln y  for the Netherlands. Figure 7 depicts the corresponding First-Order model fit (solid line) to the data (blue dots).  .
A comparison between the proposed First-Order model (19) and (20), and the Gompertz model (21), using the following linearization for the Gompertz model shows a close agreement as depicted in Figure 8. Comparisons between the proposed First-Order model and the other models, Richards, Logistic, Stannard, and Schnute, are shown in Appendix B. Table 3 summarizes the results for these comparisons by reporting the residuals obtained between the natural logarithms of the actual data and the natural logarithms of each model's fitted data. These parameters were obtained via weighted linear least squares, when linearization of the models was possible, and via nonlinear leastsquares methods, when linearization was not possible, which included modified Gauss-Newton, and Nelder-Mead (downhill simplex method) techniques and other appropriate methods for the minimization of the objective function, the sum of squares of the residuals or errors between the actual and fitted (or estimated) data. The lower the reported residuals mean value and/or the residuals RSS value, the better the fit. We find from the comparisons with other log(N)-type models that our model results are in very close agreement and that resulting curves are almost indistinguishable in the graphical rendering of these pair-wise comparisons between our model and the other models utilized, including those with larger number of parameters, as depicted in the graphs in Figure 8 and in those in Appendix B.
Other models such as in reference [10] were also studied but not compared, as we concentrated in the log(N) based ones. ear least squares, when linearization of the models was possible, and via nonlinear least squares methods, when linearization was not possible, which included modified Gauss-Newton, and Nelder-Mead (downhill simplex method) techniques and other appropriate methods for the minimization of the objective function, the sum of squares of the residual or errors between the actual and fitted (or estimated) data. The lower the reported resid uals mean value and/or the residuals RSS value, the better the fit. We find from the comparisons with other log(N)-type models that our model result are in very close agreement and that resulting curves are almost indistinguishable in the graphical rendering of these pair-wise comparisons between our model and the othe models utilized, including those with larger number of parameters, as depicted in the graphs in Figure 8 and in those in Appendix B.
Other models such as in reference [10] were also studied but not compared, as we concentrated in the log(N) based ones.

Discussion
The model presented here is meant to be a simple one that can be applied to infectiou diseases and other biological growth phenomena that follow a sigmoidal behavior. The list of models that have been developed and/or put to use in practice triggered by the

Discussion
The model presented here is meant to be a simple one that can be applied to infectious diseases and other biological growth phenomena that follow a sigmoidal behavior. The list of models that have been developed and/or put to use in practice triggered by the COVID19 pandemic is widely varied in their underlying mathematical approach, using statistical, biological-based, dynamical systems, time-series analysis, and others. The Center for Disease Control and Prevention (CDC) [22] receives forecasts from several contributors that include Universities, Research Groups, Defense, and Private Corporations (Google, Walmart, etc.) that provide both National and State predictions for a forecast interval of 4 weeks, for both new and total presumed deaths, and provides an ensemble forecast at the National and State level, based on those provided by 35 modeling groups. In [22] there is a link with information on, and links to, the various models [23]. Another link [24] provides predictions or forecasts comparisons between the model for the US and States, which are very comparable in their 4 weeks forecasts. One model the authors considered of interest was the one developed by the Discrete Dynamical Systems (DDS) group [25] based on statistics and on discrete linear dynamical systems (LDS) theory. Even though the complexity of the models described are quite varied, it is interesting to note key pros and cons factors on each, as described by the New England Journal of Medicine on their article [26] on epidemiological models: "Wrong but Useful-What COVID-19 Epidemiologic Models can and cannot tell us".

Conclusions
In this paper, we derived and compared a simple model for estimating the evolution of a disease and other biological growth processes that follow, as it occurs in many cases, a natural sigmoidal evolutionary response. Models can obviously be improved or expanded to accommodate any particular systemic bias or model variations in the data that is collected, but complexity comes at the cost of less robustness in the predictions, or higher parameter sensitivities (partial derivatives) unless a dynamical model is employed that continually updates the model parameters as new data is considered that follow latest data trends, or that incorporates a fading memory of earlier collected data. These dynamical models can be designed based on time-series analysis that includes forecasting and control theory and stochastic processes models. While the set of possibilities in creating models for estimation, prediction, or forecasting is quite large, it is always convenient to have a simple model or a set of proven simple models as a basis for building more complex ones.
well. Considerations related to the χ 2 fitting are utilized to derive the a posteriori value for σ. Therefore, initially, we assume a value σ for our weights, then fit the model parameters by minimizing χ 2 , and finally recomputing When the measurement error is not known (in our case the error for reported numbers), this method provides some error-bar to be assigned to the points. Taking the derivative of (A1) with respect to the parameters a k , we obtain equations that must hold for the minimum: For the case of fitting the data to a straight line (as is our case here), the approach just described follows: Consider the problem of fitting a set of N data points (x i , y i ) to the straight-line model We assume that the uncertainty σ i associated with each measurement y i is known, and that x i is known exactly. To measure how well the model agrees with the data, we apply the chi-square merit Function (A1), which for this case is For determining the parameters that minimize this function, we apply the usual rules from calculus: we can write the system (A6) as aS + bS x = S y , and aS x + bS xx = S xy .
Solving this (simple) system of equations (e.g., via Cramer's rule) yields ∆ = S · S xx − (S x ) 2 ; a = 1 ∆ S xx · S y − S x · S xy ; b = 1 ∆ S · S xy − S x · S y (A9) Next, we need to determine the uncertainties in the parameters a and b. Considering propagation of errors, the variance σ 2 f in the value of any function f will be For the straight line, the derivatives of a and b with respect to y i can be directly evaluated from their solutions (A9): Summing over the points we get

Appendix B
This Appendix depicts the comparisons of our proposed First Order Model to the Richards, Logistic, Stannard, and Schnute models, illustrated in Figures A1-A4, respectively, and the modified forms of the equations are described in Equations (A13)-(A16). This Appendix B complements the information provided in the Section 3.2, Comparison to other Models.  There are also modified forms of the Logistic, Gompertz, Richards, and Schnute models from reference [12,17,21], whose simpler forms were discussed in Section 3.2. Comparisons to other models are listed below for completeness. The parameter µ m is the slope of the tangent line at the inflection point of the sigmoid (at t = t i ), i.e., the maximum derivative value of the sigmoid.        Logistic-Modified Form: Gompertz-Modified Form: Richards-Modified Form: Schnute-Modified Form: