InTechOpen uses cookies to offer you the best online experience. By continuing to use our site, you agree to our Privacy Policy.

Computer and Information Science » Numerical Analysis and Scientific Computing » "Numerical Analysis - Theory and Application", book edited by Jan Awrejcewicz, ISBN 978-953-307-389-7, Published: September 9, 2011 under CC BY-NC-SA 3.0 license. © The Author(s).

Chapter 6

Numerical Simulation on Ecological Interactions in Time and Space

By Kornkanok Bunwong
DOI: 10.5772/24945

Article top


The number of influenza-associated pediatric deaths.
Figure 1. The number of influenza-associated pediatric deaths.
Examples of time scales.
Figure 2. Examples of time scales.
The time series solution of Equation (9) when α=3.4.
Figure 3. The time series solution of Equation (9) when α=3.4.
The bifurcation diagram of μ when α=3.6.
Figure 4. The bifurcation diagram of μ when α=3.6.
The time series solution of Equation (9) when α=3.6and μ=1.
Figure 5. The time series solution of Equation (9) when α=3.6and μ=1.
The time series solution of Equation (10) on time scale P1,1.
Figure 6. The time series solution of Equation (10) on time scale P1,1.
The time series solution of Equation (10) on time scale P0.1,1.
Figure 7. The time series solution of Equation (10) on time scale P0.1,1.

Figure 8.
Time evolution of the density of infected individuals b1=0, d1≠0. Parameters: qi and b0=0.3, b1=0, d0=0.2, d1=0(from top to bottom, respectively).
Figure 9. Time evolution of the density of infected individuals b1=0, d1≠0. Parameters: qi and b0=0.3, b1=0, d0=0.2, d1=0(from top to bottom, respectively).

Numerical Simulation on Ecological Interactions in Time and Space

Kornkanok Bunwong1

1. Introduction

The most formal and systematic tool to simplify a real-world phenomenon dealing with the interrelationships between organisms and their environment, including the interaction between each other is an ecological model. Of course, it is usually created for one or more purposes, for example, to gain system understanding, to forecast the future state of the system, and to develop new hypotheses. Not only that, a good model should also reach a balance between the complexities of the real-world system, which is too difficult to solve, and the need of simple formulation and valid analytic model, which can be explicitly solved. Ecological modeling has a long and rich history. So far many sophisticated ecological models have been developed in all fields of ecology such as the ecology of individuals including physiological ecology, the ecology of populations, and the study of ecosystems (classified by Hofbauer & Sigmund, 1988; Roughgarden, 1996). Ecological phenomenon can often be simplified by making some assumptions and studying it with suitable time scales and spatial interactions. It should be noted that a simple deterministic system could behave dramatically unlike when its time scale and its spatial interaction are changed. Therefore, the next coming sections aim to convince the significance of time scales and spatial interactions.

2. The classical concept of ecological modelling

The books on mathematical modeling in biology written by L. Edelstein-Keshet (1988), S. P. Ellner & J. Guckenheimer (2006), M. R. S. Kulenovic & O. Merino (2002), J. D. Murray (1993), and S. P. Otto & T. Day (2007) are the valuable teaching resources suitable for modelers who require theoretical concepts as well as desire mathematical and computational techniques to build and analyze models. In the early ecological models, their time scales were limited to either discrete or continuous. Usually, if populations of individuals have synchronized changes (i.e., reproduction, infection, recovery, migration, removal or mortality) at a regular time interval and no overlap between successive generations then the discrete time model is recommended. Otherwise, a continuous time model is preferred. Traditionally, a discrete time model is represented by using difference equation(s) whereas a continuous time model is constructed by using differential equation(s). Recently, the practical use of mathematical modeling is richly contributed. For example, L. J. S. Allen (1994), L. J. S. Allen & A. B. Burgin (2000), W. M. Getz & J. O. Lloyd-Smith (2006), and S. R.-J. Jang (2008) were interested in discrete epidemic models while K. Bunwong et al. (2009) studied nutrient removal process on a continuous time model.

According to spatial interactions, the early ecological models were based on the mass-action law, first coined in chemical reaction. Hence the system is homogenous and the spatial interaction neither exists nor plays an important role. Before the 1970s, mathematical modelers typically used ordinary differential equations, seeking equilibria and analyzing their stability (Neuhauser, 2001). Subsequently, exploring bifurcation diagrams and chaotic patterns has continuously received lots of attention, related articles were written by J. Awrejcewicz (1991), J. Awrejcewicz & C.-H. Lamarque (2003), and E. Ott (1993). In 1969, the implicit spatial model, known as metapopulation, was first introduced by R. Levins. He constructed a model to describe a population that consists of several sub-populations joined together with immigration and emigration. Levins's simple model was motivated by and applied to a pest control situation over a large region, within which local populations would fluctuate in asynchrony (Hanski & Gilpin, 1991). This new concept was closely linked with the processes of local extinctions and re-colonization. I. Hanski & M. Gilpin (1991) also provided a conceptual distinction between local, metapopulation, and geographical scales. Later on, I. Hanski (1999) applied this approach to conservation biology. However, the role of space at the individual level is still not directly mentioned.

Throughout this chapter, an SIS epidemic model, well known disease transmission model, is considered. Of course, there are lots of diseases in this world. Here we focus on the disease that does not produce immunity, for example, some STD’s, the eye disease, and the common cold. Therefore, the main situation is that the population is divided into a susceptible (S) group and an infectious (I) group. The S group is infected by the I group while the I group recovers from the disease and returns to the S group. Consequently, we set S(t) and I(t) to represent the number of S and I individuals at time t, respectively. The total population size is assumed to be constant, S(t)+I(t)=N. Moreover, (α/N)S(t)I(t) and γI(t) represent the infection rate at which the S population contracts the disease and the total number of I individuals who recover per unit time at the time t, respectively. Then, the SIS epidemic model should contain two equations. With constant population size, it can be reduced to a single equation. For continuous time scale, it is



In this case, a solution of Equation (1) always behaves non-oscillatory. For discrete time scale, the SIS epidemic model becomes


In this case, the behavior of the endemic solution can be very complicated as it can tend to an equilibrium point, to limit cycles, or show chaos. Obviously, time scale makes system behavior amazing. In the next section, an SIS model is extended with more different time scales. Finally, the qualitative structures are investigated.

3. The important of time scales

In more complex situation, flu virus can spread continuously while in season, mostly disappear for some period of time, and spread again in a new season. According to Fig. 1, the data shows periodic outbreaks of disease (the Influenza Division of Centers for Disease Control and Prevention, 2010). Obviously, the reasonable time scale for this model should be a combination of discrete and continuous scales. Thus the model either using difference equation(s) or differential equation(s) should be improper. Consequently, the calculus on time scales has been developed. S. Hilger first introduced this theory in order to unify continuous (R) and discrete (Z) analysis (Agarwal et al., 2002). Since then, the theory has been extended. Nowadays, time scales theory can be used to explain system dynamics not only for continuous and discrete times but also for other types of time such as period (P) and discrete jump with fixed length (hZ). There are plenty of publications in theoretical results. For example, J. Hoffacker & C. C. Tisdell (2005) studied on stability and instability for dynamic equations while E. Akin et al. (2001), D. R. Anderson (2009), and Y. Xu & Z. Xu (2009) studied on oscillation and nonoscillation criteria for dynamic equations and dynamic systems. However, there are few articles in numerical results (Sae-jie & Bunwong, 2009; Siming et al., 2008). In addition, there are some applications in economics (Atici et al., 2006) and epidemiology (Sae-jie et al., 2010a, 2010b; Thomas et al., 2009). In order to contribute time scales concepts and techniques understandably, some useful notations and definition are introduced as follows.


Figure 1.

The number of influenza-associated pediatric deaths.

3.1. Notations and definitions

A time scale, T, is an arbitrary nonempty closed subset of the real numbers. Forward and backward jump operators are defined by σ(t)=inf{sT:s>t} and ρ(t)=sup{sT:s<t}, respectively where inf=supT, sup=infT, and denotes the empty set. A point tT is called left-dense if t>infT and ρ(t)=t, right-dense if t<supT and σ(t)=t, left-scattered if ρ(t)<t, right-scattered if σ(t)>t, isolated if ρ(t)<t<σ(t), and dense if ρ(t)=t=σ(t). Moreover, the set Tκ is defined by T\{m} if T has a left-scattered maximum m. Otherwise, it is T. Finally, the graininess function μ:T[0,) is defined by μ(t)=σ(t)t (Bohner & Peterson, 2001).

3.2. Calculus on time scales

Traditionally, limits and continuity are key concepts for calculus development including calculus on time scales. A function f:TR is said to be rd-continuous (right dense continuous) provided f is continuous at right-dense points and left-hand limits exist and it is finite at left-dense points in T. Assume f:TR is a function and tTκ. Then the following statements are equivalent:

The (delta) derivative of f:TR at point tTκ, fΔ(t), exists.

  • For all ε>0, there is a neighborhood Uof t (i.e., U=(tδ,t+δ)T for some δ>0) such that for all sU,

  • |f(σ(t))f(s)fΔ(t)(σ(t)s)|ε|σ(t)s|

    Apart from discrete and continuous time scales, there are two more interesting time domain (as visualized in Fig. 2) to introduce. The first one is the combination of continuous and discrete time scales, denoted by the symbol T=Pl,h:=k=0[k(l+h),k(l+h)+l] where l,h>0 and kN0. For this period time scale, l is the fixed length of the continuous interval while h is the fixed length of the discrete jump. The second one is composed of points that are equally spaced in time. Suppose that the distant between two successive points is h. Therefore, the symbol of this time scale is T=hZ={hk:kZ,h>0}. Later on, forward and backward jump operators, graininess function, and (delta) derivative of function for each time scale are provided.


    Figure 2.

    Examples of time scales.

    In the case T=R, we have

    σ(t)=ρ(t)=t, μ(t)=0,


    In the case T=hZ, we have

    σ(t)=t+h,ρ(t)=th, μ(t)=h,


    In the case T=Pl,h, if L=k=0[k(l+h),k(l+h)+l), H0=k=0{k(l+h)+l}, and H1=k=1{k(l+h)} then we have





    For more details on time scales, we refer the reader to Bohner & Peterson (2001,2003)

    3.3. SIS epidemic model on time scales

    After replacing dSdt in Equation (3) with SΔ(t), we obtain the single equation of SIS epidemic model on time scales. Previously, one solution of continuous time model tends to one asymptotically stable equilibrium point while the solution of discrete time model displays various behaviors depending on parameter values. The following analyzes reveal system behavior on different time scale via two approaches. Firstly, we are interested in SIS epidemic model on T=hZ, i.e.,


    It should be pointed out that μ and h are equivalent from now on. μ is not a graininess function. Obviously, T=Z when h=1 and T=R when h tends to zero. Secondly, we change our focus to SIS epidemic model on T=Pl,h, i.e.,




    Then we explore the numerical results of previous dynamic equations individually. For T=hZ time scale, we particularly investigate the system behavior when h or μ vary. For T=Pl,h time scale, we fix h=1. Therefore, T=Z when l approaches zero.

    3.4. Numerical solution

    Since the numerical method for ordinary differential equation and difference equation are already well-known, this section contains only some additional algorithm for computing the numerical solution of dynamic system on four time scales (Sae-jie, & Bunwong, 2009). Before running the following program, the user must define the function F and G.

    INPUT: A: starting time

    X, Y: initial values

    M: number of step size (Case 2: M = l/H)

    H: step size for continuous interval

    Case 1:

    FOR j = 0, 1, 2,..., M DO

    T = A + H*j

    SAVE T, X, Y in the OUTPUT LIST



    Case 2:

    P: numbers of period

    l: length of continuous interval

    h: length of jump

    DUM3 = A

    FOR j = 0, 1, 2,..., M DO

    T = DUM3 + H*j

    SAVE T, X, Y in the OUTPUT LIST

    DUM1 = X

    DUM2 = Y



    X = DUM1 + h*F(T, DUM1, DUM2)

    Y = DUM2 + h*G(T, DUM1, DUM2)

    FOR k = 1, 2, 3,..., P DO

    DUM3 = A + k*(l+h)

    FOR j = 0, 1, 2,..., M DO

    T = DUM3 + H*j

    SAVE T, X, Y in the OUTPUT LIST

    DUM1 = X

    DUM2 = Y



    X = DUM1 + h*F(T, DUM1, DUM2)

    Y = DUM2 + h*G(T, DUM1, DUM2)


    Case 3:

    DUM3 = A

    FOR j = 0, 1, 2,..., M DO

    T = DUM3 + H*j

    SAVE T, X, Y in the OUTPUT LIST

    DUM1 = X

    DUM2 = Y

    X = DUM1 + H*F(T, DUM1, DUM2)

    Y = DUM2 + H*G(T, DUM1, DUM2)


    F := XΔ(T, X, Y)

    G := YΔ(T, X, Y)

    In the case T=hZ, W. Sae-jie et al. (2010a) always fixed the following parameter values γ=0.9 and N=100. For each value of α, μ was treated as the bifurcation parameter. When α=2, the dynamic behavior of continuous and discrete system are the same. For α=3.4, the solutions, however, appear as asymptotically stable, a period two cycle, a period four cycle when μ=0.1,0.9,1.0, respectively as shown in Fig. 3.


    Figure 3.

    The time series solution of Equation (9) when α=3.4.

    M. R. S. Kulenovic & O. Merino (2002) and L. Tien-yien & J. A. Yorke (1975) discovered that if there exists a period three cycle, then there exists chaotic behavior. For α=3.6, the bifurcation diagram (as shown in Fig. 4) exhibits a period three cycle for μ(1.0476,1.0524). Consequently, dynamic equation (9) can generate chaotic pattern for a proper value of μ. For more details, the non-oscillatory solution, the oscillating period two solution, and the chaos occur when μ=0.1, μ=0.8, and μ=1, respectively. The last behavior is illustrated in Fig. 5.

    In the case T=Pl,h, this period time scale and a continuous time scale are equivalent when μ=h=0. Thus, the result appears as a non-oscillatory solution. However, W. Sae-jie et al. (2010b) always fixed the length of discrete jump together with the following parameter values α=3.6, μ=h=1, γ=0.9, and N=100 but varied the length of continuous interval, l. For sufficiently high value of l, the result (as shown in Fig. 6) appears as a non-oscillatory solution which is similar to the result in Fig. 3. for a continuous time scale. It disappears in some intervals because of the discrete time jump. When the length of the continuous interval decreases, a period two cycle with some continuous intervals appears as visualized in Fig. 7. If the length of the continuous interval is gradually reduced and closely zero, T=Pl,h is similar to a discrete time scale.


    Figure 4.

    The bifurcation diagram of μ when α=3.6.


    Figure 5.

    The time series solution of Equation (9) when α=3.6and μ=1.

    4. The important of spatial interaction

    In contrast to a metapopulation model, space in an explicit model should be taken into account as another variable. As for non spatial models, the time scale, space, and population state for an explicit model can be either discrete or continuous. If all of them are continuous variables, a partial differential equation has been used (for example, Fisher, 1937). However, solving partial differential equations analytically is difficult. Consequently, analysis is based on computer simulation. Then we face another limitation because solving them numerically


    Figure 6.

    The time series solution of Equation (10) on time scale P1,1.


    Figure 7.

    The time series solution of Equation (10) on time scale P0.1,1.

    involves discretizing space and time. In principle, this reduces them to coupled map lattice model. On the other hand, if space is treated as discrete, a lattice model is proposed. With respect to state variables, lattice models can be divided into three subgroups, these are coupled map lattice (Morris, 1997), cellular automaton (Wolfram, 1986), and network model (Sole & Manrubia, 1996, 1997; Verdasca et al., 2005). A major problem with the lattice model is that analysis is often restricted to direct computer simulation which consumes expensive time when it is a stochastic model on a reasonably large lattice.

    Fortunately, in microscopic point of view, there is an alternative way, so-called pair approximation, that can help us not only to produce the numerical solutions of spatially explicit lattice models by using ordinarily numerical method but also to analyze them. The Japanese researchers were pioneers to apply this idea to ecological systems and have continuously produced theoretical result and applications (Harada et al., 1995; Iwasa et al., 1998; Kubo et al., 1996; Matsuda et al., 1992; Sato et al., 1994). They referred to the pair approximation as the doublet decoupling approximation and represent their system as coupled system of ordinary differential equation for the density of singletons and pairs. However, Rand’s colleagues (Bauch, 2000; Bunwong 2006; Keeling, 1995; Morris, 1997; van Baanlen & Rand, 1998) used pair approximation as a moment closure approximation where the number of pairs and the number of singletons are the only state variables. Higher order correlations are neglected. A more advantage of this method is that it keeps spatial correlations. There are varieties of applications in ecological interactions (Dieckmann et al., 2000; Ellner, 2001; Ellner et al., 1998), epidemiology (Benoit et al., 2006; de Aguiar et al., 2004 ; Dieckmann et al., 2002 ; Elliott et al., 2000 ; Joo & Lebowitz, 2004), and forest dynamics (Schlicht & Iwasa, 2007). However, the mathematical formulas are still limited. Thus, we attempt to broaden the pair approximation idea by calculating configuration averages. We call our method as a new approach and previous method as an original approach. In order to contribute our techniques understandably, some useful notations and definition are introduced as follows.

    4.1. Notations and definitions

    Under a given configuration σ¯=(σk) where k{x,e}, the following Rand’s notations are defined (Rand, 1998).

    σx,σeare the state of the site x and the edge e, respectively,σx=iis that the state of the individual x is i, σe=ijis that one end of the edge e is in state i, ei, while the other is in state j, ej, [i],[ij],[ijk] are the number of sites, edges, and triples in state i, ij, and ijk, respectively, Qx(i),Qej(i) are the number of i-state neighbors of the sites x and ei, respectively,Qx(i)σx=jis the average value of the number of i-state neighbors of a j-state site, Qej(i)σe=jkis the average value of the number of i-state neighbors of a j-state site in a jk-state edge,qiequals [i]/N where N is the total population size, andqi|jequals [ij]/Q[j] where Q is the average number of neighbors.

    In this framework, space is represented by a network of sites. Each site can either be occupied by an individual or remains as an empty site that is still available for an individual to occupy. Two sites are neighbors when they regularly interact with each other. Joining these two neighboring sites performs an edge or pair. A line is used for this interaction. Fig. 8 provides an example when the state of site xis i and the state of site y is j. Both sites are neighbors. Moreover, there is an edge e.


    Figure 8.

    4.2. SIS Master and correlation equations

    Our correlations are microcorrelations which can be measured on the scale of the interactions of individuals. After approximating higher order terms in master equations, we obtain a system of ordinary differential equation which composed of density of lower order terms, known as correlation equations. Moreover, the approximation technique is called the moment closure approximation. The differential equation for the single numbers involves pair numbers, triple numbers, etc. The differential equation for the pair numbers involves triple numbers, etc. So we get an infinite hierarchy of equations. However, we need to truncate the hierarchy at some point. For instant, pair approximation, the first order of moment closure approximation, truncates triples and higher order terms as functions of singletons and pairs only (Rand, 1998).

    Let f be a real-valued function of the state of the network at time t, which can be approximated as continuous. The equation f is derived by summing over all events in the population which affect f and the total change produced by those events is


    where dfdt=εeventsr(ε)Δfε is the rate of event r(ε) and ε is the change produced in f by event Δfε. It is called the master equation.

    For our case study, the state of each site and edge will change over time as a consequence of two major types of events – infection and recovery. Infection changes the state ε of the edge σe=SI into the state e at rate σ'e=II and recovery changes the state β of a site σx=I into the state x at rate σ'x=S. Therefore, the SIS spatial model becomes


    If d[II]dt=2σe=SIβQeS(I)2σx=IδQx(I) and δ are constant, the original approach of pair approximation is still valid. However, the real-world situation is more complicated. For example, the human-to-human transmission of Swine Flu occurs by inhalation of infectious droplets and droplet nuclei, and by direct contact, which is facilitated by air and land travel and social gatherings (Sinha, 2009). Therefore, the transmission rate and the recovery rate could vary depending on the surrounding infectious people. Consequently, we are able to assume that the infection rate and the recovery rate are β and β=b0+b1QeS(I), respectively where δ=d0d1Qx(I) are constant. Then the formulation of pair approximation is in trouble. Using the fact that

    b0,b1,d0,d1 and


    we, then, can use average forms instead of summation terms in the master equation as follows,


    4.3. Calculating configuration averages

    It should be pointed out that there are two types of average values in this new approach. Previously, all mentioned average values are called space or population average value because they are averages of the quantity over subsets of population. Now we introduce another way to calculate an average value. That is the expected value with respect to probability distribution. This seems reasonable to assume that if population size is large enough, then the configuration averages approximate probability expectations. Under some necessary and sufficient conditions, the space average and probability average are identical. K. Bunwong (2006, 2010a, 2010b) developed more formulas under multinomial distribution with parameters d[II]dt=2b0[SI]QeS(I)σe=SI+2b1[SI]QeS(I)QeS(I)σe=SI2d0[I]Qx(I)σx=I+2d1[I]Qx(I)Qx(I)σx=I and Q where pi and Poisson distribution for coordination numbers, respectively. In this chapter, the SIS spatial model is based on the framework that each site connects to a fixed number of neighbors, pi=qi|j. The following Bunwong’s formulas (2006) are used.


    4.4. Numerical solution

    The numerical method for spatial model (14) is as same as for ordinary differential equation. Here, the number of pairs and the number of singletons are state variables. K. Bunwong (2010a) mainly investigated the density of infected individuals, defined by Qej(l1)Qej(l2)σe=ij={(Q1)ql1|j+(Q1)!(Q3)!ql1|j2;l1=l2(Q1)!(Q3)!ql1|jql2|j;l1l2=qi, along the time series. We always fixed [I]/N. In case that the infection rate and the recovery rate are not affected by the surrounding infectious individuals (b0=0.3,d0=0.2), Fig. 9 reveals that the more nearby neighbors, the higher density of infected individuals at the equilibrium point. Moreover, in case that only the infection rate is affected by the surrounding infectious individuals (b1=0,d1=0), the stronger effect of the surrounding infectious individuals on the infection rate, the higher density of infected individuals at the equilibrium point. In case that only the recovery rate is affected by the surrounding infectious individuals (b10,d1=0), the stronger effect of the surrounding infectious individuals on the recovery rate, the higher density of infected individuals at the equilibrium point. The illustrations can be seen in Bunwong (2010a).


    Figure 9.

    Time evolution of the density of infected individuals b1=0,d10. Parameters: qi and b0=0.3,b1=0,d0=0.2,d1=0(from top to bottom, respectively).

    5. Conclusion

    Obviously, the assumption on timescale affects the system behavior. In order to simplify real-world correctly, choosing a suitable time scale is important for model formulation. Our example, SIS epidemic model, proves that if continuous time scale is used then two solutions of the system are asymptotically stable or unstable depending on parameter values and stable oscillating solutions have never existed. In contrast, if discrete time scale is applied then there are various types of solution behaviors such as equilibrium point solutions, period two cycles, period four cycles, period three cycles, and also chaotic solutions depending on parameter values as well. Consequently, the predicted behaviors from a model can be qualitatively very different for different time scales. The hard to answer question is “What is the proper model to understand observed data in a variety of time measurements?”. Of course, the observed data is usually discrete. Can we use differential equation(s)? With theory of time scales, the models are more varieties. Bifurcation diagram shows that the discrete jump and the continuous interval are essential. The differential equation(s) can be described observed data if the jump distant of data is sufficiently small. The distant also depends on other parameters. Moreover, numerical verification can visualize other interesting behavior pattern and guide other mathematician to consider theoretically.

    In recent years, the effect of spatial structure is often taken into account in ecological interactions. Pair approximation is one of the powerful tools to understand human to human interactions. We have developed the way to approximate higher order quantities and applied to ecological problems. Particularly, our new approach is suitable for a model evolving according to the transition rates affecting additionally by neighbors. For example, people infect flu virus easily from their nearby neighbors. The health organization usually suggests infectious people to have some rest and be away from public places. It implies that if we surrounding with more infectious people, then we have higher chance to infected and/or lower chance to recover as shown in the numerical results.

    Further studies on time scales should involve more numerical techniques and applications while future works on spatial interaction should concern development of formulation.


    I would like to thank Prof. David Rand for introducing me to spatial world, Prof. Yoh Iwasa for sharing his experience, Prof. Ravi Agarwal for introducing me to time scale world, Dr. Elvin Moore for intuitive ideas, and Dr. Wichuta Sae-jie for collaboration. I am grateful to the Development and Promotion of Science and Technology Talent Project (DPST) and Thailand Research Fund (TRF) for financial supports under contact number MRG5180246.


    1 - R. Agarwal, M. Bohner, D. O’Regan, A. Peterson, 2002Dynamic Equations on Time Scales: A Survey. Journal of Computational and Applied Mathematics, 1411-21260377-0427
    2 - E. Akin, L. Erbe, A. Peterson, B. Kaymakçalan, 2001Oscillation Results for a Dynamic Equation on a Time ScaleJournal of Difference Equations and Applications767938101023-6198
    3 - L. J. S. Allen, 1994Some Discrete-Time SI, SIR, and SIS Epidemic Models. Mathematical Biosciences, 1241831050025-5564
    4 - L. J. S. Allen, A. B. Burgin, 2000Comparison of Deterministic and Stochastic SIS and SIR Models in Discrete Time. Mathematical Biosciences, 16311330025-5564
    5 - D. R. Anderson, 2009Oscillation and Nonoscillation Criteria for Two-Dimensional Time-Scale Systems of First-Order Nonlinear Dynamic Equations. Electronic Journal of Differential Equations, 241131072-6691
    6 - F. M. Atici, D. C. Biles, A. Lebedinsky, 2006An Application of Time Scales to Economics. Mathematical and Computer Modelling, 437-87187260895-7177
    7 - J. Awrejcewicz, 1991Bifurcation and Chaos in Coupled OscillatorsWorld Scientific Publishing, 9-81020-579-1
    8 - J. Awrejcewicz, C. H. Lamarque, 2003Bifurcation and Chaos in Nonsmooth Mechanical SystemsWorld Scientific Publishing, 9-81238-459-6
    9 - C. T. Bauch, 2000Moment Closure Approximations in Epidemiology.Warwick PhD Thesis.
    10 - J. Benoit, A. Nunes, M. T. da, Gama, 2006Pair Approximation Models for Disease Spread. European Physical Journal B, 501-21771811434-6028
    11 - M. Bohner, A. Peterson, 2001Dynamic Equations on Time Scales: An Introduction with Applications,Birkhäuser, 0-81764-225-0USA.
    12 - M. Bohner, A. Peterson, 2003Advances in Dynamic Equations on Time ScalesBirkhäuser, 0-81764-293-5USA.
    13 - K. Bunwong, 2006Spatial Modeling in Evolutionary Ecology.Warwick PhD Thesis.
    14 - K. Bunwong, W. Sae-jie, Y. Lenbury, 2009Modelling Nitrogen Dynamics of a Constructed Wetland: Nutrient Removal Process with Variable Yield. Nonlinear Analysis: Theory, Methods & Applications, 7112e1538e15460036-2546X.
    15 - K. Bunwong, 2010Pair Approximations for Ecological Models with Additional Neighborhood EffectsProceedings of the 4th International Conference on Applied MathematicsSimulation, Modelling (ASM ‘10), 165169978-9-60474-210-3Corfu Island, Greece.
    16 - K. Bunwong, 2010A New Approach to Pair Approximation Method for Spatial Model in Ecology. WSEAS Transactions on Mathematics, 9107687771109-2769
    17 - M. A. M. de Aguiar, E. M. Rauch, Y. Bar-Yam, 2004Invasion and Extinction in the Mean Field Approximation for a Spatial Host-Pathogen Model. Journal of Statistical Physics, 1145-6141714510022-4715
    18 - U. Dieckmann, R. Law, J. A. J. Metz, 2000The Geometry of Ecological Interactions: Simplifying Spatial Complexity,Cambridge University Press, 0-52164-294-9UK.
    19 - U. Dieckmann, J. A. J. Metz, M. W. Sabelis, K. Sigmund, 2002Adaptive Dynamics of Infectious Diseases: In Pursuit of Virulence Management,Cambridge University Press, 0-52178-165-5UK.
    20 - L. Edelstein-Keshet, 1988Mathematical Models in BiologyRandom House, 0-39435-507-5York, USA.
    21 - P. Elliott, J.C. Wakefield, N.G. Best, D.J. Briggs, (2000) Spatial Epidemiology: Methods and Applications, Oxford University Press, ISBN 0-1926-2941-7, New York, USA
    22 - S. P. Ellner, 2001Pair Approximation for Lattice Models with Multiple Interaction Scales. Journal of theoretical biology, 21044354470022-5193
    23 - S.P. Ellner, J. Guckenheimer, (2006). Dynamic Models in Biology, Princeton University Press, ISBN 0-6911-1843-4, New Jersey, USA
    24 - S. P. Ellner, A. Sasaki, Y. Haraguchi, H. Matsuda, 1998Speed of Invasion in Lattice Population Models: Pair-Edge ApproximationJournal of Mathematical Biology3654694840303-6812
    25 - R. A. Fisher, 1937The Wave of Advance of Advantageous Genes. Annals of Human Genetics. 74355369
    26 - W. M. Getz, J. O. Lloyd-Smith, 2006Basic Methods for Modeling the Invasion and Spread of Infectious Diseases. DIMACS Series in Discrete Mathematics and Theoretical Computer Science, 711231052-1798
    27 - Y. Harada, H. Ezoe, Y. Iwasa, H. Matsuda, K. Sato, 1995Population Persistence and Spatially Limited Social Interaction. Theoretical population Biology, 48165910040-5809
    28 - I. Hanski, (1999). Metapopulation Ecology, Oxford University Press, ISBN 0-1985-4065-5, New York, USA. (Reprinted 2005
    29 - I. Hanski, M. Gilpin, 1991Metapopulation Dynamics: Brief History and Conceptual Domain. Biological Journal of the Linnean Society, 421-23161095-8312Reprinted in Gilpin, M. & Hanski, I. (eds). Metapopulation Dynamics: Empirical and Theoretical Investigations, Academic Press, 0-12284-120-4
    30 - J. Hofbauer, K. Sigmund, 1988The Theory of Evolution and Dynamical Systems, Cambridge University Press,0-52135-838-8UK.
    31 - J. Hoffacker, C. C. Tisdell, 2005Stability and Instability for Dynamic Equations on Time Scales. Computers and Mathematics with Applications, 499-10132713340898-1221
    32 - Y. Iwasa, M. Nakamaru, S. A. Levin, 1998Allelopathy of Bacteria in a Lattice Population: Competition between Colici-Sensitive and Colicin-Producing Strains. Evolutionary Ecology, 127858020269-7653
    33 - S. R. J. Jang, 2008Backward Bifurcation in a Discrete SIS Model with Vaccination. Journal of Biological Systems, 1644794940218-3390
    34 - J. Joo, J. L. Lebowitz, 2004Pair Approximation of the Stochastic Susceptible-Infected-Recovered-Susceptible Epidemic Model on the Hypercubic Lattice. Physial Review E., 7031539-3755
    35 - M. J. Keeling, 1995The Ecology and Evolution of Spatial Host-Parasite SystemsWarwick PhD Thesis.
    36 - T. Kubo, Y. Iwasa, N. Furumoto, 1996Forest Spatial Dynamics with Gap Expansion: Total Gap Area and Gap Size Distribution. Journal of Theoretical Biology, 18032292460022-5193
    37 - M. R. S. Kulenovic, O. Merino, 2002Discrete Dynamical Systems and Difference Equations with MathematicaChapman and Hall/CRC, 1-58488-287-5USA.
    38 - H. Matsuda, N. Ogita, A. Sasaki, K. Sato, 1992Statistical Mechanics of Population-the Lattice Lotka-Voltera Model. Progress of Theoretial Physics, 88103510490003-3068X.
    39 - A. J. Morris, 1997Representing Spatial Interactions in Simple Ecological Models.Warwick PhD Thesis.
    40 - J. D. Murray, 1993Mathematical BiologySpringer-Verlag, 354057204Berlin, German.
    41 - C. Neuhauser, 2001Mathematical Challenges in Spatial Ecology. Notices of the AMS. 481113041314
    42 - E. Ott, 1993Chaos in Dynamical SystemsCambridge University Press, 0-52143-799-7UK.
    43 - S.P. Otto, T. Day, (2007). A Biologist’s Guide to Mathematical Modeling in Ecology and Evolution. Princeton University Press, ISBN 978-0-6911-2344-8, New Jersey, USA
    44 - D. A. Rand, 1998Correlation Equations and Pair Approximations for Spatial Ecologies. Warwick Preprint.
    45 - J. Roughgarden, (1996). Theory of Population Genetics and Evolutionary Ecology: an Introducation. Prentice Hall, ISBN 0-1344-1965-0, New Jersey, USA
    46 - W. Sae-jie, K. Bunwong, 2009Numerical Verification of Certain Oscillation Results on Time Scales. Neural, Parallel & Scientific Computations, 173173381061-5369
    47 - W. Sae-jie, K. Bunwong, E.J. Moore, (2010). Qualitative Behavior of SIS Epidemic Model on Time Scales. Proceedings of the 4th International Conference on Applied Mathematics, Simulation, Modelling (ASM ‘10), pp. 159-164, ISBN 978-960-474-210-3. Corfu Island, Greece, July 22-25, 2010
    48 - W. Sae-jie, K. Bunwong, E. J. Moore, 2010The Effect of Time Scales on SIS Epidemic Model. WSEAS Transactions on Mathematics, 9107577671109-2769
    49 - K. Sato, H. Matsuda, A. Sasaki, 1994Pathogen Invasion and Host Extinction in Lattice Structured Populations.Journal of Mathematical Biology3232512680303-6812
    50 - R. Schlicht, Y. Iwasa, 2007Spatial Pattern Analysis in Forest Dynamics: Deviation from Power Law and Direction of Regeneration WavesEcology Research, 2221972030912-3814
    51 - Z. Siming, L. Yanhui, R. Zhixin, 2008Region Qualitative Analysis of Predator-Prey System on Time Scales.Annals of Differential Equations241211261002-0942
    52 - M. Sinha, (2009). Swine Flu. Journal of Infection and Public Health, Vol. 2, pp. 157-166, ISSN 1876-0341
    53 - R.V. Sole, S.C. Manrubia, (1996). Extinction and Self-Organized Criticality in a Model of Large-scale Evolution. Physical Review E, 0106-3651154R42R45
    54 - R. V. Sole, S. C. Manrubia, 1997Criticality and Unpredictability in Macroevolution. Physical Review E, 55445004507
    55 - Y. Takenaka, H. Matsuda, Y. Iwasa, 1997Competition and Evolutionary Stability of Plants in a Spatially Structured HabitatResearches on Population Ecology39167750034-5466
    56 - The Influenza Division of Centers for Disease Control and Prevention. (May 2010) Fluview: A Weekly Influenza Surveillance Report, In: Seasonal Influenza (Flu), 07.04.2011, Available from:
    57 - D. Thomas, M. Weedermann, L. Billings, J. Hoffacker, R. A. Washington-Allen, 2009When to Spray: a Time-Scale Calculus Approach to Controlling the Impact of West Nile VirusEcology and Society142Article 21, 1708-3087
    58 - L. Tien-yien, J. A. Yorke, 1975Period Three Implies Chaos. The American Mathematical Monthly, 82109659920002-9890
    59 - M. van Baanlen, D.A. Rand, (1998). The Unit of Selection in Viscous Populations and the Evolution of Altruism. Journal of Theoretical Biology, 0022-51934193631648
    60 - J. Verdasca, M.M.T. Da Gama, A. Nunes, N. R. Bernardino, J. M. Pacheco, M. C. Gomes, 2005Recurrent Epidemics in Small World Networks. Journal of Theoretical Biology, 23345535610022-5193
    61 - S. Wolfram, 1986Theory and Applications of Cellular Automata. World Scientific Publishing. 9-97150-123-6
    62 - Y. Xu, Z. Xu, 2009Oscillation Criteria for Two-Dimensional Dynamic Systems on Time Scales. Journal of Computational and Applied Mathematics, 22519190377-0427