On the Relationship of Steady State of Continuous and Discrete Models From

Dynamic Programming

Martin L. Puterman , in Encyclopedia of Physical Science and Technology (Third Edition), 2003

IV.C.3 Continous-Time Models

Stochastic continuous time models are categorized according to whether the state space is continuous or discrete. The discrete time model has been widely studied in the operations research literature. The stochastic nature of the problem is modeled as either a Markov process, a semi Markov process, or a general jump process. The decision maker can control the transition rates, transition probabilities, or both. The infinite horizon versions of the Markov and semi Markov decision models are analyzed in a similar fashion to the discrete-time Markov decision process; however, general jump processes are considerably more complex. These models have been widely applied to problems in queuing and inventory control.

When the state space is continuous and Markovian assumptions are made, diffusion processes are used to model the transitions. The decision maker can control the drift of the system or can cause instantaneous state transitions or jumps. The discrete-time optimality equation is replaced by a nonlinear second order partial differential equation and is usually solved numerically. These models are studied in the stochastic control theory literature and have been applied to inventory control, finance, and statistical modeling.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B0122274105001873

Control Systems, Identification

R.K. pearson , in Encyclopedia of Physical Science and Technology (Third Edition), 2003

I.B.1 Continuous vs. Discrete Time

Physical systems typically evolve more or less smoothly in continuous time, and fundamental models that attempt mechanistic descriptions of these systems generally consist of collections of ordinary differential equations, partial differential equations, integral equations, and so forth. Conversely, computer control is generally implemented at discrete time instants, based on measurements made at discrete time instants. Consequently, computer-based control is commonly based on discrete-time dynamic models, possibly obtained as approximate discretizations of continuous-time models. Also, a point that is particularly relavent to this discussion is that discrete-time models are generally much easier to identify from input–output data than continuous-time models. In the case of simple linear models, an exact correspondence between continuous-time and discrete-time models is possible, but this correspondence does not extend to the nonlinear case.

As an important example, consider the first-order linear dynamic model in continuous-time, described by the following ordinary differential equation:

(4) $\frac{\text{<mi mathvariant="italic" is="true">dy</mi>}(\text{t})}{\text{<mi mathvariant="italic" is="true">dt</mi>}}=-\frac{1}{\tau }\text{y}(\text{t})+\frac{\text{G}}{\tau }\text{u}(\text{t})$

Here, u(t) represents an input variable, y(t) is the corresponding output variable, τ is the time constant for the system, and G is the steady-state gain of this system. If the variables u(t) and y(t) are sampled uniformly in time and we write u _k   = u(t _k ) and y _k   = y(t _k ) where t _k   = t ₀  + kT, it is possible to obtain the following difference equation relating u _k and y _k , also known as a discrete-time dynamic model:

(5) ${\text{y}}_{\text{k}}=\alpha {\text{y}}_{\text{k}-1}+\beta {\text{u}}_{\text{k}-1}$

The basis on which Eq. (5) is derived from Eq. (4) is the assumption that u(t) is piecewise constant and does not change between sampling instants, a reasonable assumption in cases where the sequence {u _k } is computer generated and a reasonable approximation in many other cases. Under this assumption, the coefficients of the continuous-time and discrete-time models are related by:

(6) $\alpha ={\text{e}}^{-\text{T}/\tau },\beta =\text{G}(1-{\text{e}}^{-\text{T}/\tau })$

Similar equivalences exist between higher order linear ordinary differential equations and linear difference equations (Pearson, 1999b, p. 12), but it is also important to note that the transformation just described from the continuous-time model, Eq. (4), to the discrete-time model, Eq. (5), imposes significant restrictions on the parameter α. In particular, for a stable first-order linear model, the response time τ is necessarily positive; one consequence is that responses to step changes in u(t) from u ₀ for t  ≤   0 to u ₊ for t  >   0 result in responses that settle out essentially to their asymptotic value Gu ₊ in ∼3 to ∼5 time constants τ. In contrast, if τ   <   0   in Eq. (4), the step response is unstable, growing exponentially without bound as t  →   ∞. This stability range (0   <   τ   <   ∞) corresponds to the range 0   <   α   <   1 for the parameter α in the discrete-time model, whereas unstable models with τ   <   0 correspond to α   >   1   in the discrete-time model. It is important to note that the discrete-time model also exhibits a stable response for −1   <   α   <   0, but this response is oscillatory, corresponding to behavior that is only possible for linear continuous-time models of order 2 or higher. This result illustrates that continuous-time and discrete-time models can behave somewhat differently even for simple linear models where these differences are, relatively speaking, quite small.

For nonlinear models, the differences between continuous-time and discrete-time behavior are much more pronounced. This point may be seen by considering the following simple nonlinear reactor model, discussed in more detail by Pearson (1999b, Sec. 8.5):

(7) $\frac{\text{<mi mathvariant="italic" is="true">dy</mi>}(\text{t})}{\text{<mi mathvariant="italic" is="true">dt</mi>}}=-\text{h}{\text{y}}^2(\text{t})-\frac{\text{y}(\text{t})\text{u}(\text{t})}{\text{V}}+\frac{\delta \text{u}(\text{t})}{\text{V}}$

where h, V, and δ are constants. Assuming as in the previous example that u(t) is piecewise constant and defining μ _k   = u(t _k )/2hV ultimately lead to the following exact discretization:

$\begin{array}{cc}\hfill {\text{y}}_{\text{k}}& =\frac{[1-{\tau }_{\text{k}-1}{\mu }_{\text{k}-1}]{\text{y}}_{\text{k}-1}+2\delta {\tau }_{\text{k}-1}{\mu }_{\text{k}-1}}{1+{\tau }_{\text{k}-1}[{\text{y}}_{\text{k}-1}+{\mu }_{\text{k}-1}]}\hfill \\ \hfill {\tau }_{\text{k}-1}& =\frac{\text{tanh}\left[\text{<mi mathvariant="italic" is="true">hT</mi>}\sqrt{{\mu }_{\text{k}-1}^2+2\delta {\mu }_{\text{k}-1}}\right]}{\sqrt{{\mu }_{\text{k}-1}^2+2\delta {\mu }_{\text{k}-1}}}\hfill \end{array}$

The key point here is that, although the order remains the same on discretization (i.e., a first-order nonlinear differential equation results in a first-order nonlinear difference equation), the forms of these two equations are nothing alike; the simple polynomial appearing on the right-hand side of the differential equation is mapped into a messy combination of hyperbolic functions, square roots, and fractions.

It is, of course, possible to approximate a set of nonlinear differential equations by a structurally similar set of nonlinear difference equations. As a specific example, the Euler approximation simply replaces derivatives with finite differences; dy(t)/dt ≃ [y _k   − y _k−1]/T, but this approximation is not without its difficulties. As a specific example, consider the Euler discretization of the second-order nonlinear differential equation:

(8) $\frac{\text{<mi mathvariant="italic" is="true">dx</mi>}(\text{t})}{\text{<mi mathvariant="italic" is="true">dt</mi>}}=-\alpha {\text{x}}^2(\text{t})\to {\text{x}}_{\text{k}}={\text{x}}_{\text{k}-1}-\alpha \text{T}{\text{x}}_{\text{k}-1}^2.$

Comparing the responses of these equations to the initial condition x(t)   = x ₀ at t  =   0 illustrates an important difference: the exact solution of the continuous-time equation is x(t)   = x ₀/(1   +   α x ₀ t), which is stable for all x ₀  >   0 and all α   >   0, decaying asymptotically to zero like 1/αt as t  →   ∞. In contrst, the discretized model can be shown to be unstable if αx ₀ T  >   1. A comparison of the continuous-time solution and the approximate discrete-time solution is shown in Fig. 1 for α x ₀ T  =   1.026.

FIGURE 1. Continuous- vs. discrete-time model responses.

Because they are necessary for computer-based control applications and because the range of techniques available for discrete-time model identification appears to be wider than that for continuous-time model identification, discrete-time system identification methods appear to be more popular in practice and are discussed more extensively in the literature. Consequently, the remainder of this article restricts consideration to discrete-time model identification. For a general overview and some useful references on the continuous-time case, consult the survey paper by Unbehauen and Rao in Sawaragi and Sagara (1997, pp. 973–999); this survey deals with both linear and nonlinear continuous-time dynamic models, and other papers in the same conference proceedings provide useful illustrations of these and related ideas, especially those presented in the session on continuous-time identification (Sawaragi and Sagara, 1997, pp. 1281–1310). In addition, the discussion given here is also restricted to time-invariant or shift-invariant linear models, again because these models are most widely used in practice. The topic of time-varying linear model identification is also treated in a dedicated session at the SYSID'97 meeting, beginning with a useful overview by Forssell and Ljung (Sawaragi and Sagara, 1997, pp. 1655–1678).

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B012227410500140X

Special Volume: Computational Methods for the Atmosphere and the Oceans

Jacques Blum , ... I. Michael Navon , in Handbook of Numerical Analysis, 2009

5.1 Optimal nudging specification

We assume that the model equations have been discretized in space by a finite difference, finite element, or spectral discretization method. The time continuous model satisfies dynamical equations of the form

(5.1) $\frac{\partial X}{\partial t}=F\left(X\right),$

(5.2) $X\left(0\right)=V,$

where X represents the discretized state variable of the model atmosphere, t is time, and V represents the initial condition for the model. Say, for instance, X^o (t) is a given observation, then the objective of VDA is to find model initial conditions that minimize a cost function defined by

(5.3) $J\left(V\right)={\int }_0^T〈W\left(X-X^{\circ }\right),X-X^{\circ }〉\text{d}t,$

where W is a diagonal weighting matrix. Note that J is only a function of the initial state because X is uniquely defined by the model equations (Eqs. (5.1) and (5.2)). An implicit assumption made in VDA is that the model exactly represents the state of the atmosphere. However, this assumption is not true.

The NDA technique introduced by Anthes [1974] consists in achieving a compromise between the model and the observations by considering the state of the atmosphere to be defined by

(5.4) $\frac{\partial X}{\partial t}=F\left(X\right)+G\left(X^{\circ }-X\right),$

where G is a diagonal matrix.

Together with the initial conditions

(5.5) $X\left(0\right)=V,$

the system (Eq. (5.1)) has a unique solution X(V, G). The main difficulty in the NDA scheme resides in the estimation of the nudging coefficient G( Stauffer and Seaman [1990]). If G is too large, the fictitious diffusion term will completely dominate the time tendency and will have an effect similar to replacing the model data by the observations at each time-step. Should a particular observation have a large error that prevents obtaining a dynamic balance, an exact fit to the observation is not required since it may lead to a false amplification of observational errors. On the other hand, if G is too small, the observation will have little effect on the solution. In general, G decreases with increasing observation error, increasing horizontal and vertical distance separation, and increasing time separation. In the experiment of Anthes [1974], a nudging coefficient of 10⁻³was used for all the fields for a hurricane model and was applied on all the domain of integration. In the experiment of Krishnamurti, Jishan, Bedi, Ingles and Oosterhof [1991], the relaxation coefficients for the estimated NDA experiment were kept invariant both in space and time, and their values were simply determined by numerical experience. The implicit dynamic constraints of the model then spread the updated information to the other variables (temperature and moisture) resulting eventually in a set of balanced conditions at the end of the nudging period.

In the work of Zou, Navon and Le Dimet [1992], a new parameter estimation approach was designed to obtain optimal nudging coefficients. They were optimal in the sense that the difference between the model solution and the observations will be small. For a comprehensive review of parameter estimation addressing issues of identifiability, see Navon [1998].

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S1570865908002093

Simulation of Dynamic Models

Mark M. Meerschaert , in Mathematical Modeling (Fourth Edition), 2013

6.2 Continuous-Time Models

In this section we discuss the fundamentals of simulating continuous-time dynamical systems. The methods presented here are simple and usually effective. The basic idea is to use the approximation

$\frac{dx}{dt}\approx \frac{\Delta x}{\Delta t}$

to replace our continuous-time model (differential equations) by a discrete– time model (difference equations). Then we can use the simulation methods we introduced in the preceding section.

Example 6.2. Reconsider the whale problem of Example 4.2. We know now that starting at the current population levels of B = 5, 000, F = 70, 000, and assuming a competition coefficient of α < 1.25 × 10⁻⁷, both populations of whales will eventually grow back to their natural levels in the absence of any further harvesting. How long will this take?

We will use the five-step method. Step 1 is the same as before (see Fig. 4.3), except that now the objective is to determine how long it takes to get to the equilibrium starting from B = 5, 000, F = 70, 000.

Step 2 is to select the modeling approach. We have an analysis question that seems to require a quantitative method. The graphical methods of Chapter 4 tell us what will happen, but not how long it will take. The analytical methods reviewed in Chapter 5 are local in nature. We need a global method here. The best thing would be to solve the differential equations, but we don't know how. We will use a simulation; this seems to be the only choice we have.

There is some question as to whether we want to adopt a discrete-time or a continuous-time model. Let us consider, more generally, the case of a dynamic model in n variables, x = (x ₁,…, x_n ), where we are given the rates of change F = (f ₁,…, f_n ) for each of the variables x ₁,…, x_n , but we have not yet decided whether to model the system in discrete-time or continuous-time. The discrete-time model looks like

(6.4) $\begin{array}{l}\Delta x_1=f_1\left(x_1,\dots ,x_n\right)\hfill \\ ⋮\hfill \\ \Delta x_n=f_n\left(x_1,\dots ,x_n\right),\hfill \end{array}$

where Δx_i represents the change in x_i over 1 unit of time (Δt = 1). The units of time are already specified. The method for simulating such a system was discussed in the previous section.

If we decided on a continuous-time model instead, we would have

(6.5) $\begin{array}{l}\frac{{dx}_1}{dt}=f_1\left(x_1,\dots ,x_n\right)\hfill \\ ⋮\hfill \\ \frac{{dx}_n}{dt}=f_n\left(x_1,\dots ,x_n\right),\hfill \end{array}$

which we would still need to figure out how to simulate. We certainly can't expect the computer to calculate x(t) for every value of t. That would take an infinite amount of time to get nowhere. Instead we must calculate x(t) at a finite number of points in time. In other words, we must replace the continuous– time model by a discrete-time model in order to simulate it. What would the discrete-time approximation to this continuous-time model look like? If we use a time step of Δt = 1 unit, it will be exactly the same as the discrete-time model we could have chosen in the first place. Hence, unless there is something wrong with choosing Δt = 1, we don't have to choose between discrete and continuous. Then we are done with step 2.

Step 3 is to formulate the model. As in Chapter 4, we let x ₁ = B and x ₂ = F represent the population levels of each species. The dynamical system equations are

(6.6) $\begin{array}{ll}\frac{{dx}_1}{dt}\hfill & =0.05x_1\left(1-\frac{x_1}{150,000}\right)-\alpha x_1{x}_2\hfill \\ \frac{{dx}_2}{dt}\hfill & =0.08x_2\left(1-\frac{x_2}{400,000}\right)-\alpha x_1{x}_2\hfill \end{array}$

on the state space x ₁ ≥ 0, x ₂ ≥ 0. In order to simulate this model we will begin by transforming to a set of difference equations

(6.7) $\begin{array}{l}\Delta x_1=0.05x_1\left(1-\frac{x_1}{150,000}\right)-\alpha x_1{x}_2\hfill \\ \Delta x_2=0.08x_2\left(1-\frac{x_2}{400,000}\right)-\alpha x_1{x}_2\hfill \end{array}$

over the same state space. Here, Δx_i represents the change in population x_i over a period of Δt = 1 year. We will have to supply a value for α in order to run the program. We will assume that α = 10⁻⁷ to start with. Later on, we will do a sensitivity analysis on α.

Step 4 is to solve the problem by simulating the system in Eq. (6.7) using a computer implementation of the algorithm in Fig. 6.2. We began by simulating N = 20 years, starting with

$\begin{array}{ll}{x}_1\left(0\right)\hfill & =5,000\hfill \\ x_2\left(0\right)\hfill & =70,000.\hfill \end{array}$

Figures 6.9 and 6.10 show the results of our first model run. Both blue whale and fin whale populations grow steadily, but in 20 years they do not get close to the equilibrium values

Figure 6.9. Graph of blue whales x ₁ versus time n for the whale problem: case α = 10⁻⁷, N = 20.

Figure 6.10. Graph of fin whales x ₂ versus time n for the whale problem: case α = 10⁻⁷, N = 20.

$\begin{array}{ll}{x}_1\hfill & =35,294\hfill \\ x_2\hfill & =382,352\hfill \end{array}$

predicted by our analysis back in Chapter 4.

Figures 6.11 and 6.12 show our simulation results when we input a value of N large enough to allow this discrete-time dynamical system to approach equilibrium.

Figure 6.11. Graph of blue whales x ₁ versus time n for the whale problem: case α = 10⁻⁷, N = 800.

Figure 6.12. Graph of fin whales x ₂ versus time n for the whale problem: case α = 10⁻⁷, N = 100.

Step 5 is to put our conclusions into plain English. It takes a long time for the whale populations to grow back: about 100 years for the fin whale, and several centuries for the more severely depleted blue whale.

We will now discuss the sensitivity of our results to the parameter α, which measures the intensity of competition between the two species. Figures 6.13 through 6.18 show the results of our simulation runs for several values of α. Of course, the equilibrium levels of both species change along with α.

Figure 6.13. Graph of blue whales x ₁ versus time n for the whale problem: case α = 3 × 10⁻⁸, N = 800.

Figure 6.14. Graph of blue whales x ₁ versus time n for the whale problem: case α = 10⁻⁸, N = 800.

Figure 6.15. Graph of blue whales x ₁ versus time n for the whale problem: case α = 10⁻⁹, N = 800.

Figure 6.16. Graph of fin whales x ₂ versus time n for the whale problem: case α = 3 × 10⁻⁸, N = 100.

Figure 6.17. Graph of fin whales x ₂ versus time n for the whale problem: case α = 10⁻⁸, N = 100.

Figure 6.18. Graph of fin whales x ₂ versus time n for the whale problem: case α = 10⁻⁹, N = 100.

However, the time it takes our model to converge to equilibrium changes very little. Our general conclusion is valid whatever the extent of competition: It will take centuries for the whales to grow back.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780123869128500063

Modeling Formalisms and Their Simulators

Bernard P. Zeigler , ... Ernesto Kofman , in Theory of Modeling and Simulation (Third Edition), 2019

3.2.2 Continuous System Simulation

After describing the solutions of LTI systems, we are back to the problem of finding the solution of a general ODE:

(3.18) $\dot{\mathbf{z}}(t)=\mathbf{f}(\mathbf{z}(t),\mathbf{u}(t))$

with known initial state $\mathbf{z}(t_0)={\mathbf{z}}_0$ .

As we already mentioned, in most cases the exact solution cannot be obtained so we are forced to obtain numerical approximations for the values of the state at certain instants of time $t_0,t_1,\cdots ,t_N$ . These approximations are obtained using numerical integration algorithms, and traditionally the term Continuous System Simulation is tightly linked to those algorithms.

Nowadays, Continuous System Simulation is a wider topic. In fact, most continuous time models are not originally written like Eq. (3.18). That representation is the one used by numerical integration routines but, for a modeling practitioner, it is not comfortable (and it is sometimes impossible) to express a model in that way.

Modern continuous systems modeling languages like Modelica – the most widely accepted modeling standard language – allow representing the models in an object oriented fashion, reusing and connecting components from libraries. The resulting model is then an object oriented description that, after a flattening stage, results in a large set of differential and algebraic equations. Then, those Differential Algebraic Equations (DAEs) must be converted in an ODE like that of Eq. (3.18) in order to use numerical integration routines. This conversion requires different (and sometimes complex) algorithms and it constitutes itself an important branch of the Continuous System Simulation field (an extensive discussion is in Chapter 3 of TMS2000).

Many continuous time models describe the evolution of some magnitudes both in time and space. These models are expressed by Partial Differential Equations (PDEs) and their simulation constitute one of the most difficult problems in the discipline.

In this book we shall not cover the problems of DAE or PDE simulation. However, as both DAEs and PDEs can be converted or approximated by ODEs, the numerical algorithms for ordinary differential equations described next can still be applied for those problems.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128133705000110

Transmission of Infectious Diseases

Winfried Just , ... Natalia Toporikova , in Algebraic and Discrete Mathematical Methods for Modern Biology, 2015

8.2.5 How to Model Time and Run Simulations

Running simulations is somewhat similar to observing many outbreaks in a tightly controlled artificial environment. Making sense of simulation output is a little bit like drawing conclusions from real-world data. In particular, the observed results will share some common features in different runs of the simulation, but there also may be substantial variability from run to run due to the inherent stochasticity of disease transmission. In Section 8.2 of [17], you will learn how to use our software and how to deal with this inherent messiness so that the common features can be clearly discerned.

However, simulation outputs are also very much unlike real data. The environment in which they run is tightly controlled both by the underlying model that is embodied in the code and by the input parameters. It is important to investigate how much the simplifying assumptions that go into the process of building the model might influence its predictions, as reflected in the simulation outputs. In Chapter 9, we will investigate in detail how assumptions about the mixing pattern influence the predictions of the resulting models.

When investigating a given model, one still needs to choose its parameters, such as the mean duration 〈τ ^I 〉 of infectiousness. Epidemiologists need to estimate these parameters from real data, but may not be able to do so with great precision. Moreover, control measures often translate into a change of some of these parameters. It is therefore important to study how the outcomes of simulations change in response to changes in the parameters. Section 8.3 of [17] gives you an introduction to this type of investigation, which goes by the name of sensitivity analysis. It contains several relevant simulation exercises. You can work through these exercises at any time after completing Section 8.2 of [17].

A related issue is the fact that the population size in simulations is usually rather limited (to a few hundred hosts in our code), but real populations can have millions of hosts. It is therefore important to study how the predictions of simulations scale up, that is, how they change when the population size N becomes large. This is a special case of sensitivity analysis (sensitivity to the parameter N), but the goal here is different: epidemiologists usually know the size of the target population with good precision and are not trying to design control measures that will modify it, ⁶ but are interested in how likely simulations performed on a small population would give reliable information about outbreaks in large populations. We will return to this issue in Chapter 9.

It is time to run some actual simulations. We strongly recommend that before reading the remainder of this chapter you work through the exercises of Section 8.2 of [17].

When working through the simulation exercises of Section 8.2 of [17], you may have noticed the mysterious drop-down menu model-time of IONTW. We perceive time as advancing continuously, and the setting model-time $\to$ Continuous makes the computer treat time in this way. In terms of our models, this means that state transitions could be happening at any time, as they will in reality. For example, the time $T_{13}^I$ at which host 13 becomes infectious could take any one of the values 3,4,3.1,3.2,3.14,3.15,…. We could even have $T_{13}^I=\pi$ .

But suppose the onset of infectiousness of host 13 "really" occurs at time $T_{13}^I=\pi$ and time is measured in days. We could presumably observe that the onset of infectiousness occurred during the fourth day, that is, $3<T_{13}^I=\pi \le 4$ , or perhaps even that it happened between 2:24 a.m. and 4:48 a.m. On the fourth day, which translates into $3.1<T_{13}^I\le 3.2$ . But could we possibly tell whether $3.14<T_{13}^I\le 3.15$ ? Note that the difference between 3.14 and 3.15 days is 14.4 min. It is becoming clear that we will never be able to obtain data that pinpoint the transition times with arbitrary precision. An equation like $T_{13}^I=\pi$ is only a mathematical abstraction, albeit a very useful one.

Exercise 8.12

Our definition of transition events allows for the possibility that, for example, $T_6^I=T_{77}^I=T_8^R=t$ for some t. In other words, each of these state transitions would happen exactly at the same time. Why can we ignore this possibility in continuous-time models?

In continuous-time models, the disease transmission parameters need to be specified in terms of two rates β and α. If an E-compartment is included in the model, we also need a third rate constant γ, but to keep matters as simple as possible we will disregard this option in the discussion of this section. The parameter β represents the rate at which two hosts make effective contact, and it allows us, for any given Δt ≥ 0, to calculate the probabilities b _i,j that we introduced in Section 8.2.1. The parameter α represents the rate at which hosts cease to be infectious. It allows us, for any given Δt ≥ 0, to calculate the probability a that an infectious host ceases to be infectious during an interval of length Δt as well as the mean duration of infectiousness 〈τ ^I 〉. Details of these calculations will be explained in Section 8.6 of [17].

Often it is easier to pretend that time progresses in discrete, integer-valued steps, and to simply consider the state transitions that are going to happen until the next time step. Note that the epidemic curves of Figure 8.1 of Section 8.1 implicitly treat time in this way: each time step in these graphs corresponds to 1 day, and the curves supposedly tell us about the number of new infections that occurred during a given day. For day t + 1 this would be the total number of hosts i for which $t<T_i^I\le t+1$ . For example, for day 3 one would set t = 2 and count that total number of hosts i for which $2<T_i^I\le 3$ .

We can go one step further and base our simulations on discrete-time models. In these models, time is supposed to take only integer values t = 0,1,2,…, and instead of considering exact transition times $T_i^I$ and $T_i^R$ , we work with the rounded-up versions $\lceil T_i^I\rceil$ and $\lceil T_i^R\rceil$ . Recall that the ceiling function ⌈T⌉ is defined so that if t < T ≤ t + 1, then ⌈T⌉ = t + 1. Time t = 0 will be used only for designating the initial state. At subsequent times t + 1 we proceed as follows. If $\lceil T_i^I\rceil =t+1$ , then we move host i into the I-compartment. If $\lceil T_i^R\rceil =t+1$ , then we move host i into the R-compartment (in the case of an SIR-model) or back into the S-compartment (in the case of an SIS-model).

Exercise 8.13

Consider the state transition sequence (8.2). Find the state transition events for times t = 1,2,3 in the corresponding discrete-time model.

In these models, we lose some information, though. For example, $\lceil T_i^R\rceil =t+1$ tells us only that $t<T_i^R\le t+1$ ; we no longer know whether removal of infectious host i occurred near the beginning or near the end of the time interval (t,t + 1]. This might matter though, as in the latter case host i will still be able to infect other hosts during the time interval (t,t + 1], while in the former case this is very unlikely. We will explore shortly what this loss of information does to the predictions of our model.

When constructing discrete-time models, we need to carefully decide how to choose the length of the interval of physical time that corresponds to one time step of the discrete model. If we consider an interval of 1 week, then t would signify the number of weeks that have elapsed since the start of the outbreak. This may be a good choice if, for example, the mean duration of infectiousness 〈τ ^I 〉 is 1 month, so that hosts will stay infectious for about 4 consecutive time steps. But if, on average, hosts stay infectious only for 2 days, this would be a bad choice, as most hosts would become infectious and recover during the same time interval, that is, $t<T_i^I<T_i^R\le t+1$ for some t. In an SIS-model where hosts move back into the S-compartment upon recovery, we would not even see most infections.

If the option model-time $\to$ Discrete is chosen, the input field time-step of IONTW allows you to explore various choices of the time scale for discrete-time modeling. We will explore the effects of these choices in Exercise 8.14 below.

What if the mean duration of infectiousness 〈τ ^I 〉 is exactly 1 week? There is always some variability, so some hosts will still both become infectious and recover during the same time interval, and for some hosts there will be a time t with $T_i^I<t<t+1<T_i^R$ . But if the variability is small, these two scenarios may happen so infrequently that we might be able to ignore them. At time t + 1 we could then move all hosts who reside in the I-compartment at time t into the R-compartment (in the case of an SIR-model) or back into the S-compartment (in the case of an SIS-model). We will refer to such greatly simplified discrete-time models as next-generation models. In Section 8.5 of [17], you can find a detailed explanation where this terminology comes from. For simulating such models with our software IONTW, the parameter end-infection-prob that represents the probability a of removal by the next time step needs to be set to 1.

It might appear that these models are useful for predicting the time course of an outbreak only if the mean duration of infectiousness happens to be a natural time unit like 1 week or 1 month. But they could work just fine if the mean duration of infectiousness is, for example, 10.4 days. All we would need to do is rescale our interpretation of physical time to units of 10.4 days; unfamiliar to be sure, but not essentially different from rescaling days into weeks. Our next exercise will illustrate to what extent we might be able to trust the predictions of next-generation and other discrete-time models.

Exercise 8.14

Compare the predictions of four SIR-models of the same disease. In one model, choose the continuous-time option; in the second and third models, choose the discrete-time approximation with time steps Δt = 0.02 and Δt = 1, respectively. In the fourth model, set the probability of removal a to 1 so as to obtain a next-generation model. Detailed instructions for the exercise are given in Section 8.4 of [17].

Take-Home Message:

•

Continuous-time models keep track of the exact moments at which state transition events occur. Time t is represented by nonnegative real numbers. In simulations of continuous-time models without an E-compartment, the parameters are specified in terms of two rates β and α. In our software IONTW, the input fields infection-rate and end-of-infection-rate control the values of these parameters.

•

In discrete-time models, time t is represented by nonnegative integers that reflect scaling of time in suitable units Δt. All transition events that occur in the interval (t,t + 1] of real time are treated as if they occurred simultaneously at time t + 1. In simulations of discrete-time models, the parameters are specified in terms of transition probabilities and the length of the time step Δt. In our software, the input fields infection-prob and end-of-infection-prob control the values of these parameters, while time-step controls the value of the parameter Δt.

•

Next-generation models are discrete-time models based on time units Δt = 〈τ ^I 〉. In these models, it is assumed that infectious hosts cease to be infectious after exactly one time unit. In order to simulate such models with our software, we need to set the value of end-of-infection-prob to 1 and the value of time-step to 〈τ ^I 〉.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128012130000083

Digital Controllers

David M. Auslander , in Encyclopedia of Physical Science and Technology (Third Edition), 2003

VI.D Difference Equations

The state space and transfer function formats describe the continuous-time portion of the system, normally the control object. In computer-controlled systems, however, the controller is discrete-time and the control object is continuous-time. The discrete-time portion of the system is described by difference equations rather than the differential equations that describe the continuous-time portion. The discrete-time portion of the system description has two portions: a part that has been converted from the continuous-time model and a part for the controller. The first part is converted from a linear model, so will also be linear. The second part, the controller, is specified by the system designer. It is "linearized" by omitting the nonlinear parts (such as saturations).

The state space format for the discrete time equations is very similar in appearance to the continuous-time version,

(7) $\begin{array}{lll}\hfill \mathbf{x}(\text{k}+1)& =\hfill & \mathbf{Fx}(\text{k})+\mathbf{Gu}(\text{k})\hfill \\ \hfill \mathbf{y}(\text{k})& =\hfill & \mathbf{Cx}+\mathbf{Du}\hfill \end{array}$

The coefficient matrices are sometimes given different names to avoid confusion, and sometimes given the same names to indicate similarity of function. The index k is the sample time counter.

Because the controller is normally the "open" part of the design, it is customary to convert the entire system to the discrete-time domain. There are a number of methods for converting the continuous-time part to discrete-time. If the equations can be solved (or the transforms inverted) the conversion can be done exactly ("exactly" means that the responses will match at the sampling instants; nothing is said about the behavior between samples). There are also several approximations that are simpler to use, but only valid for specific ranges of sampling intervals.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B0122274105001757

Mathematical Modeling

Xavier J.R. Avula , in Encyclopedia of Physical Science and Technology (Third Edition), 2003

III Classification of Mathematical Models

Classification of models is important to the understanding of a system's behavioral data and to interpreting or communicating the data for prediction and decision making. An understanding of the system behavior is a contribution to science, and its use for prediction and decision making constitutes technology. Mathematical models are broadly classified into dynamic and steady-state models. Obviously, the element of time is present in the dynamic models. The dynamic models are further classified into continuous-time and discrete-time models. In a continuous-time model the time advances smoothly through real numbers, while in the discrete-time model the time advances in finite jumps, each jump representing a multiple of a selected time unit.

Models constructed in terms of descriptive state variables are classified into continuous state, discrete-state, and mixed-state models. In a continuous-state model the range sets of state variables are presented by real numbers, while in a discrete-state model the variables assume a discrete set of values. In a mixed-state model both kinds of variables are present.

The broad class of continuous-time models is further subdivided into models in which the state changes continuously and those in which the state changes in discrete jumps. The former subclass of models are represented by ordinary and/or partial differential equations, for example, for the boundary-layer flow over a wing of an aircraft. On the other hand, the production line in a factory is a discrete-state system.

Yet another classification breaks down mathematical models into stochastic (or probabilistic) and deterministic models. Models with at least one random variable in their description are termed stochastic. Models other than stochastic are deterministic. In fact, deterministic is a special case of stochastic.

Another classification of models is based on how the real system interacts with the environment. If the real system is isolated from its environment, the model is called autonomous. If the environment exerts any kind of influence through so-called input variables from the environment without being controlled by the model, then the model is called nonautonomous.

When models are expressed in terms of mathematical equations, the solutions are profoundly affected by nonlinearity. Based on the type of equations, the models are classified as linear or nonlinear. These classifications can be combined with the above-mentioned model types to yield, for example, linear or nonlinear dynamic and linear or nonlinear stochastic models. Although the universe is predominantly nonlinear, some processes are linear within certain bounds. Actually, linearity can be viewed as a special case of nonlinearity. Sometimes, it is mathematically expedient to consider nonlinear processes as piece-wise linear. Systems are classified as linear if and only if they simultaneously satisfy the principle of homogeneity and the principle of superposition. The principle of homogeneity preserves the input function scale factor in transition from input to output. A system satisfies the superposition principle when the superposition of individual input functions result in a response that is the superposition of corresponding outputs. Simulations and exact solutions of linear systems of equations are well established and reported in mathematical literature.

All systems that are not linear belong to the class of nonlinear systems. As there are no general methods of solutions to solve and analyze models of nonlinear systems, simulation has become a commonly used tool for this purpose. Impacted by recent developments in numerical analysis and by advances in powerful, high-speed and high-memory digital and analog computers, simulation of nonlinear systems have blossomed into technical endeavors that are broad in scope and variety encompassing biological, ecological, social, and economic systems in addition to those in hard scientific fields such as physics, chemistry, material science, and mechanics. Modern nonlinear equations frequently encountered in simulation include, for example, the quaternion rotational equations of motion that are treated with geometric algebra which is of late applied to a range of problems in many fields of science and being promoted as a unified mathematical language in the 21st century.

Traditional quantitative techniques of modeling cannot be effectively applied to model phenomena in the so-called "soft" disciplines consisting of humanistic systems such as social, economic, ecological, and biological systems because there may arise difficulties associated with multidimensionality, subsystem interactions, inexplicable feedback mechanisms, hierarchical structures, and unpredictable behavioral dynamics. Such difficulties may also be encountered in mechanistic systems that may include some physical and engineering systems. These difficulties lead to making imprecise statements, introducing fuzziness much like human reasoning, about the characteristics and behavior of a system. An important class of models using fuzzy sets was launched about 35 years ago. The idea of fuzzy sets is centered on the imprecision in the belonging of elements to a set. Intuitively, a fuzzy set is a set of elements that are not precise and in which there is no distinct boundary between the elements that belong to the set and those that do not. In other words, the transition from full membership to nonmembership of an element (or elements) is blurred, so that a gradation of partial membership is possible. Mathematical concepts based on this idea have been successfully applied to modeling of systems in "soft" desciplines as well as in engineering when imprecision is present. Fuzzy systems models are categorized into two types: liguistic and rule-based. An elaborate discussion of these types of models and modeling techniques is beyond the scope of this article.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B0122274105004117

Nonparametric Methods in Continuous-Time Finance: A Selective Review *

Zongwu Cai , Yongmiao Hong , in Recent Advances and Trends in Nonparametric Statistics, 2003

3.1 Nonparametric estimation of parametric diffusion models

As is well-known, derivative pricing in mathematical finance is conveniently tractable in continuous-time. In empirical study, however, it is an usual practice to abandon continuous-time modeling when estimating derivative pricing models. This is mainly because the transition density for most continuous-time models with discrete observations has no closed form and so the maximum likelihood estimation (MLE) is infeasible.

One major focus of the continuous-time literature is on developing methods to estimate continuous-time models using discretely-sampled data. This is motivated by the fact that using the discrete version of a continuous-time model can result in inconsistent parameter estimates ([71]). Available estimation procedures include the MLE method of [71], the simulated methods of moments of [40] and [50], the generalized method of moments (GMM) of [53], the efficient method of moments (EMM) of [48], the Markov chain Monte Carlo (MCMC) of [61], and the methods based on the empirical characteristic functions of [63,82].

Below we focus on some nonparametric. estimation methods of a parametric model

(10) $d{X}_t=\mu \left(X_t,\theta \right)dt+\sigma \left(X_t,\theta \right)d{B}_t\text{,}$

where μ(·,·) and σ(·,·) are of known form, and θ is an unknown parameter vector in a parameter space Θ. [4] proposes a minimum distance estimator

(11) $\widehat{\theta }=arg\underset{\theta \in \mathrm{\Theta }}{min}{n}^{-1}{\displaystyle \sum _{\mathit{\tau }=1}^n{\left[{\widehat{\pi }}_0\left(X_{\mathit{\tau }\mathrm{\Delta }}\right)-\pi \left(X_{\mathit{\tau }\mathrm{\Delta }},\theta \right)\right]}^2\text{,}}$

where ${\widehat{\pi }}_0\left(x\right)=n^{-1}{\mathrm{\Sigma }}_{\mathit{\tau }=1}^n{K}_h\left(x-{X_{\mathit{\tau }}}_{\mathrm{\Delta }}\right)$ is a kernel estimator for the marginal density of X_t , and

(12) $\pi \left(x,\theta \right)=\frac{c\left(\theta \right)}{{\sigma }^2\left(x,\theta \right)}exp\left[{\displaystyle {\int }_{x_0^{·}}^x\frac{2\mu \left(u,\theta \right)}{{\sigma }^2\left(u,\theta \right)}du}\right]\text{,}$

is the marginal density estimator implied by model (10). The factor c(θ) ensures that π(·, θ) integrates to 1 for every θ ∈ ⊖, and x*₀ is the lower bound of the support of {X_t }. Because the marginal density cannot capture the full dynamics of {X_t }, $\widehat{\theta }$ is not be asymptotically most efficient., although it is root-n consistent for θ ₀.

Recently, [6] proposes an approximated likelihood approach. Let p_x (Δ, x | x ₀, θ) be the conditional density of Χ _τΔ  = x given Χ _{(τ  −   1)Δ}  = x ₀ induced by model (10). The log-likelihood function of the model for the sample {Χ _τΔ} ⁿ _{τ  =   1} is

$l_n\left(\theta \right)={\displaystyle \sum _{\mathit{\tau }=1}^{n}ln{p}_x\left(\mathrm{\Delta },X_{\mathit{\tau }\mathrm{\Delta }}|X_{\left(\mathit{\tau }-1\right)\mathrm{\Delta }},\theta \right)\text{.}}$

The MLE estimator that maximizes l_n (θ) is asymptotically most efficient. Unfortunately, except for some simple models, p_x (Δ, x|x ₀, θ) has no closed form. [6] uses a Hermite polynomial approximation p ^(J) _x (Δ, x|x ₀, θ) for p_x (Δ, x|x ₀, θ), and then obtains an estimator that maximizes the approximated model likelihood. This estimator enjoys the same asymptotic efficiency as the (infeasible) MLE as J  = J_n   →   ∞. More specifically, [6] considers a transform

$Z_{\mathrm{t}}={\mathrm{\Delta }}^{-1/2}\left(Y_t-y_0\right),\mathrm{where}{Y}_t\equiv \gamma \left(X_t,\theta \right)={\displaystyle {\int }_{-\infty }^{X_{\mathrm{t}}}{\sigma }^{-1}\left(u,\theta \right)du\text{,}}$

and then approximates the transition density of Z_t by the Hermite polynomials:

$p_z^{\left(J\right)}\left(\mathrm{\Delta },z|z_0,\theta \right)=\varphi \left(z\right){\displaystyle \sum _{j=0}^J{\eta }_z^{\left(j\right)}\left(z_0,\theta \right)H_j\left(z\right)\text{,}}$

where ϕ(·) is the N(0, 1) density, and {H_j (z)} is the Hermite polynomial series. The coefficients {η ^(j) _z (z ₀, θ)} are specific conditional moments of process Z_t , and can be calculated using the Monte Carlo method or a higher Taylor series expansion in Δ. The approximated transition density of X_t is then given as follows:

$p_x\left(\mathrm{\Delta },x|x_0,\theta \right)={\mathrm{\Delta }}^{-1/2}{p}_z\left\{{\mathrm{\Delta }}^{-1/2}\left[\gamma \left(x,\theta \right)-\gamma \left(x_0,\theta \right)\right]|\gamma \left(x_0,\theta \right),\theta \right\}\text{.}$

Under suitable conditions, the estimator ${\widehat{\theta }}^{\left(J\right)}=arg{min}_{\theta \in \mathrm{\Theta }}{\mathrm{\Sigma }}_{t=1}^{n}ln{p}_x^{\left(J\right)}\left(\mathrm{\Delta },X_{\mathit{\tau }\mathrm{\Delta }}|X_{\left(\mathit{\tau }-1\right)\mathrm{\Delta }},\theta \right)$ is asymptotically equivalent to the infeasible MLE. [5] applies this method to estimate a variety of diffusion models for spot interest rates, and finds that J  =   2 or 3 gives accurate approximations for most financial diffusion models. [42] extend this approach to stationary time-inhomogeneous diffusion models, [7] to general multivariate diffusion models and [8] to affine multi-factor term structure models.

In a rather general continuous-time setup which allows for stationary multi-factor diffusion models with partially observable state variables (e.g., stochastic volatility model), [48] propose an EMM estimator that also enjoys the asymptotic efficiency as the MLE. The basic idea of EMM is to first use a Hermite-polynomial based semi-nonparametric (SNP) density estimator to approximate the transition density of the observed state variables. This is called the auxiliary SNP model and its score is called the score generator, which has expectation zero under the model-implied distribution when the parametric model is correctly specified. Then, given a parameter setting for the multi-factor model, one may use simulation to evaluate the expectation of the score under the stationary density of the model and compute a chi-square criterion function. A nonlinear optimizer is used to find the parameter values that minimize the proposed criterion.

Specifically, suppose {X_t } is a stationary possibly vector-valued process with the conditional density p ₀(Δ, Χ _τΔ| X _sΔ, s  ≤ τ  −   1)   = p ₀(Δ, Χ _τΔ| Y _τΔ), where Y _τΔ  =   (Χ _{(τ  −   1)Δ},…, Χ _{(τ  −   1)Δ})′ for some fixed integer d  ≥   0. This is a Markovian process of order d. To estimate parameters in model (10) or its multivariate extension, [48] propose to check whether the following moment condition holds:

(13) $M\left({\beta }_n,\theta \right)\equiv {\displaystyle \int \frac{\partial logf\left(\mathrm{\Delta },x,y;{\beta }_n\right)}{\partial {\beta }_n}}p\left(\mathrm{\Delta },x,y;\theta \right)\mathit{dxdy}=0,\mathrm{if}\theta ={\theta }_0\in \mathrm{\Theta }\text{,}$

where p(Δ, x, y;θ) is the model-implied joint density for (X _τΔ, Y ′ _τΔ))′ θ ₀ is the unknown true parameter value, and f(Δ, x, y;β_n ) is an auxiliary SNP model for the joint density of (X _τΔ, Y′ _τΔ)′ Note that β_n is the parameter vector in f(Δ, x, y;β_n ) and may not nest parameter θ. By allowing the dimension of β_n to grow with the sample size n, the SNP density f(Δ, x, y;β_n ) will eventually span the true density p ₀(Δ, x, y) of (X _τΔ, Y ′ _τΔ)′, and thus is free of misspecification asymptotically. [48] use a Hermite polynomial approximation for f(Δ, x, y;β_n ), with the dimension of β_n determined by such model selection criteria as BIC. Because p(Δ, x, y;θ) usually has no closed form, the integration in (13) can be computed by simulating a large number of realizations under model (10).

The EMM estimator is defined as follows:

$\widehat{\theta }=arg\underset{\theta \in \mathrm{\Theta }}{min}M{\left({\widehat{\beta }}_n,\theta \right)}^{\prime }{\widehat{I}}^{-1}\left(\theta \right)M\left({\widehat{\beta }}_n,\theta \right)\text{.}$

where $\widehat{\beta }$ is the quasi-MLE for β_n , the coefficients in the SNP density model f(x, y;β_n ) and the matrix $\widehat{I}\left(\theta \right)$ is an estimate of the asymptotic variance of $\sqrt{n}\partial M_n\left({\widehat{\beta }}_n,\theta \right)/\partial \theta$ (see [49]). This estimator $\widehat{\theta }$ is asymptotically as efficient as the (infeasible) MLE.

The EMM has been applied widely in financial applications. See, for example, [1,15,37] for interest rate applications, [13,32,69] for estimating stochastic volatility models for stock prices with such complications as long memory and jumps, [33] for estimating and testing target zero models of exchange rates, and [64] for price option pricing. It would be interesting to compare the EMM method and approximate MLE of [6] in finite samples.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444513786500193

MAXIMISATION PRINCIPLES IN EVOLUTIONARY BIOLOGY

A.W.F. Edwards , in Philosophy of Biology, 2007

5 EVOLUTION AS FITNESS-MAXIMISATION. (3) WHAT REMAINS TO BE SAID

If it is true that the course of evolution cannot in general be described in terms of a potential function like those that exist in physics, is it possible to salvage anything from the idea at all in the context of evolution?

First, a distinction must be made between complete 'dynamic' theories involving a potential function that can fully describe the processes of change, and weaker 'static' theories of maximisation which merely describe stationary points of a system in terms of the maximisation of some non-decreasing function. The latter class includes the basic multiallelic single-locus model of population genetics in which stable equilibria are at a maximum of the mean fitness, as proved by a number of authors around 1960 (see [Edwards, 2000a], for the history). However, since the proof is specific to the problem no general theory emerges.

Somewhere between the dynamic and the static theories there is the possibility that changes can be described by a variational principle, but not one that involves a potential function. This turns out to be the case for the above model provided a transformed gene-frequency space is adopted, as follows.

Svirezhev [1972] stated that in the continuous-time model for selection at a multiallelic locus the application of Fisher's angular transformation ${sin}^{-1}\sqrt{p}$ to the gene frequencies results in (1) the direction of the gene-frequency change being given by the steepest direction on the function W and (2) the magnitude of the change being proportional to the square root of the additive genetic variance, all with respect to the new coordinate system. This does not in itself define a potential function, though the case of only two alleles exhibits an interesting atypical feature in the discrete-generation case because the square root of the additive genetic variance is exactly proportional to the slope (taken positively) of the mean fitness W in the transformed space [Edwards, 2000b]. Thus the magnitude of the change is indeed proportional to the slope and W is then a potential function.

It is possible to pursue this line of approach using the Riemannian distance Σ[(dp_i )²/p_i ] implied by Fisher's transformation, but it requires advanced mathematics and the reader is referred to the books by Svirezhev and Passekov [1990], Hofbauer and Sigmund [1998] and Bürger [2000]. However, in a more accessible account Ewens [1992] placed the result in the context of an optimization principle originally proposed by Kimura [1958] and developed by Ewens to prove that natural selection changes the gene frequencies so as to maximise the partial increase in the mean fitness subject to a constraint on the Riemannian distance between new and old gene frequencies. This is indeed a maximisation principle, but does not go so far as to establish a potential function.

We thus see that a naive description of evolution as a process that tends to increase fitness is misleading in general, and hill-climbing metaphors are too crude to encompass the complexities of Mendelian segregation and other biological phenomena. In some particular cases there is a deep structure involving variational principles, but it is principally of mathematical interest only. Fisher's Fundamental Theorem of Natural Selection seems to be about as much as it is possible to say about the relationship between genetical variation and evolutionary change, but it does not lead to any grand theory such as is to be found in the physical sciences.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444515438500174

allenwitiontis45.blogspot.com

Source: https://www.sciencedirect.com/topics/mathematics/continuous-time-model

Share this post