# A Stochastic Demographic Model as an Autonomous Non-Markov Queueing System with an Unlimited Number of Servers

by

(2017). A Stochastic Demographic Model as an Autonomous Non-Markov Queueing System with an Unlimited Number of Servers. In Young Scientist USA, Vol. 10 (p. 76). Auburn, WA: Lulu Press.

*Abstract**. In this article a mathematical model of human population growth** **as an autonomous non-Markov queuing system with an unlimited number of servers is proposed. The research is represented with a new virtual phase method and with a method of an asymptotic analysis of a stochastic age-specific density for a number of customers served in the system at time t. A probability distribution for the stochastic density of customers and main characteristics of this distribution are found. The mathematical model and methods of its research are applied to analysis of the population growth in the United State of America.*

Keywords: *mathematical model, population, queueing system, population growth, projection.*

**1. Introduction**

Demographic estimates and projections are an important element in the management of social processes. At present it is very difficult to find an area of economics or a social life where their long-term results are not used. A mathematical modeling is one of the useful methods to estimate a long-term population size.

The development and utilization of different kinds of mathematical models serve both of the population are the linear growth model and the exponential growth model [1-2], the stable population model [3-6] and its modifications [7-18], the logistic growth model [1-3, 6, 19], the hyperbolic growth model [20-21] and the Markovian models in the form of Markov chains [22-23], the cohort-component method [24-26].

The analysis of many models has shown that in the modeling of demographic processes the deterministic discrete models, the deterministic continuous models and the stochastic discrete (simulation) models are the most developed. Since the real processes of human population growth are stochastic and continuous. Thereby, the construction of a new stochastic demographic model with a continuous time is actual. Our work is devoted to find a solution to this problem. In this article we apply the models and the methods of queuing theory to estimate the long-term population size [27-28].

**2. The autonomous non-Markov queuing system with an unlimited number of servers**

We propose a mathematical model of a human population growth, for example, women in the United State of America, as an autonomous non-Markov queuing system with an unlimited number of servers (Figure 1).

**Figure 1.** The autonomous non-Markov queuing system with an unlimited number of servers

We define the process of functioning of this queuing system in terms used in queuing theory. After that we will give a demographic interpretation of the system’s main elements.

Each customer from incoming stream occupies a free server. Service durations of various customers are stochastically independent and have the same distribution, defined by a function *S(x)=1-B(x)*, where *B(x) *is the distribution function of the customers’ service duration. After finishing the service, the customer leaves the queuing system.

We define the age *x≥0* for the customer in the system as the length of the interval from the moment *t-x* of its service start (the moment of service start) till the current time *t*. Each customer at age *x* at time *t* with intensity *b(x,t)* generates new customers. Then, *b(x,t)**Δ**t+o(**Δ**t) *is the probability the customer at age *x* from the moment *t* for an infinitesimal time interval *Δ**t* generates a new customer. The probability of generating two or more customers is an infinitesimal quantity of a higher order than *Δ**t*. At the moment of appearance, a new customer occupies a free service.

In terms used in demography, a served customer is interpreted as a woman, the service duration of the customer is the life span of this woman. The function *S(x)* is the survival function for women [6], the age *x* of the customer is the age of a woman at the instant of time. The function *b(x,t)* – the fertility function of girls. The incoming flow of customers is the process of girls’ birth [6].

In the proposed stochastic continuous time model, we will investigate an age-specific stochastic density of a female population and a probability distribution of its values.

Let’s denote N(x_{1},x_{2},t) – the number of customers with the age , served in this system at time *t*. Here *x* is any nonnegative real number. A limit

will be called the stochastic density *N*(*x*,*t*) of the number of customers (female population) at age *x* at time *t*.

Obviously,

In this work we will solve the problem of the investigating the random function *N*(*x*,*t*). In demography the mathematical expectation of the random function *N*(*x*,*t*), that is

.

It is called a population density function at age *x* at time *t *[6]. The function *m*(*x*,*t*) defines the average characteristics and it is the most important characteristic of the random function N(x,t). However, in this work other characteristics of the function N(x,t), such as a variance of N(x,t), a cross-correlation function of N(x1,t) and N(x2,t)*, *a probability distribution of N(x,t) for a fixed value of *x *can be also found.

3. A research of the autonomous non-Markov queuing system with an unlimited number of servers

We perform the investigation of the random function *N(x,t)* according to a new virtual phase method. That is, we considerate an auxiliary autonomous queuing system with a PH-distribution of the service time (a phase-type distribution).

3.1. The auxiliary autonomous queuing system with the PH-distribution of the service time and its research

Consider a queuing system with a structure as in Figure 1, but the service duration τ of each customer consists of the durations of the random number of phases:

where τ* _{i}* – the duration of the

*i*

^{th}service phase. The values τ

*are independent exponentially distributed random variables with the parameter μ, where i=1,2…v. Here v is a random variable and v=1,2…. The service of each new customer begins in the first phase. The customer, finished the service in the*

_{i}*i*

^{th}phase, with probability r

*passes to the service on i+1*

_{i}^{th}phase, and with probability 1-r

*completes its service and leaves the system. Such a system will be called the autonomous system with the phase distribution or the PH-distribution of the service time.*

_{i}The number of phases v is determined by the achievement of the Markov chain (given by the transition probabilities graph) from the first step of the zero-step (Figure 2).

**Figure 2.** The transition probabilities graph of the Markov chain

We determine the transition probability of customer to the next phase as

, (1)

we believe that S(0)=1-B(0)=1. So we can write by the multiplication theorem on probability [29]

.

According to this expression and (1) it follows

, (2)

where P(v=n) is the probability that the customer service will be completed in n phase in the system. Let’s find a characteristic function g(α)=Me^{-}^{ατ} for the random variable τ. We write by the formula of total probability for mathematical expectations and according to (2)

.

For μ →∞, we rewrite this equation in the form

.

Thus, for μ→∞, the duration of the service with the PH-distribution converges to the service time determined by the distribution function *B*(*x*)=1-*S*(*x*) of the original autonomous queuing system.

We will assume that each customer served at the *i*^{th} phase at time t with intensity b* _{i}*(t)=b(i/μ,t) generates a new customer. We denote n(i,t) – the number of customer, served at the

*i*

^{th}phase at time t . Then a random process

is multidimensional the continuous time Markov chain. Its probability distribution is

.

We write the system of Kolmogorov differential equations by the formula of total probability [29]

(3)

We denote a characteristic function of the number of occupied servers at time t in the form

, * * (4)

where * – *the imaginary unit*.*

We multiply (3) on and sum over all n* _{i}*. We obtain a equality

, (5)

where , .

The solution H(u,t) of (5) determines the problem of researching the autonomous system with the PH-distribution of the service time.

**3.2. A basic equation for the autonomous non-Markov queuing system with an unlimited number of devices **

We denote

, , (6)

and the characteristic function H(u,t) of the form (4) is written as

Let a expression i/μ→x for μ→∞, i→∞, then we assume that the following limits exists

, (7)

, (8)

and there is also a distribution limit

. (9)

The function F(u,t) in the form (8) is called the characteristic functional of the random function ξ(x,t) of two arguments x and t. The random function ξ(x,t) is called the stochastic density of the number of customers of age *x* served in the system at time *t*, for all *x*≥0.

The values of the functional F(u,t) are determined by the form of the function u=u(x) and by the values of the scalar argument *t*.

Now let’s make replacement (1) in (5). For μ→∞, i/μ→x, we rewrite equation (5) by (6-9) in the form

. (10)

The equation (10) is called the basic equation for the autonomous non-Markov queuing system with an unlimited number of devices.

We will solve the equation (10) will be solved by the method of asymptotic analysis [30], assuming that the random function ξ(*x*,*t*) takes sufficiently large values proportional to some infinitely large value of *N, *that is *N* →∞. We find an asymptotic of the first and second orders of its solution for equation (10).

**3.3. A solution of the basic equation by the method of asymptotic analysis**

In queuing theory, the method of asymptotic analysis is the study of equations, determined any characteristics of a system, when an asymptotic condition is satisfied. The form of an asymptotic condition is specified in accordance with the model and the research problem [30].

We consider equation (10) for the characteristic function *H*(*u*,*t*) in the asymptotic condition of a large number of groups that are proportional to an infinitely large value of *N*. In the asymptotic condition of a large number of groups, we find the characteristic function of the number of customers *n*(*t*), serviced at time *t* [30].

**3.3.1. The first order asymptotic**

Denote ε=1/N, where ε – a small positive parameter [31], in the equation (10) we replace

, , (11)

and we obtain an equality

. (12)

This equality (12) is singularly perturbed equation [31].

We consider a subclass of solutions *H(*u,t) of equation (10), for which, by (11), the function F_{1}(w,t,ε) has the following properties.

There is a finite limit at ε→0 for F_{1}(w,t,ε)

, (13)

and its derivative with respect t and w(x)

, (14)

. (15)

We prove the following theorem.

**Theorem 1.** For ε→0 the solution of problem (12) has the form

, (16)

where m(x,t) is a solution of equation

(17)

with the boundary condition

.

**Proof.**** **In equality

**(12) the exponential is expanded in series up to**

*О*(ε

^{2}). In view of (13-15), we make the passage to the limit for ε→0. We obtain that the functional

*F*

_{1}(

*w*,

*t*) is the solution of a linear homogeneous equation

. (18)

The solution of this equation is found in the form of the characteristic functional

, (19)

where function m(x,t) is the mean stochastic age-specific density ξ(x,t) of the asymptotic number of the customer, serviced in the system at the moment t, – the imaginary unit. Substitute the solution (19) in equation (18)

.

Equating terms for the same w(0) and w(x), we obtain the equation that coincides with (17) and the corresponding boundary condition.

By virtue of (11) and also of (16), we can write the asymptotic equality for ε→0

.

**Definition***.*** **A function

is called the first order asymptotic of the characteristic functional F(u,t).

We note that if we use the asymptotic results for a demographic analysis to choose 1 year as an unit of time and *N* as 1,000, then the resulting numerical values will be expressed in thousands, as it is usually done in demographic studies. And in terms of demography, equation (17) together with the boundary condition determines the function Nm(x,t). Nm(x,t) is an average value of the stochastic population density.

The analytic solution of the partial differential equation (17) is founded by the composition and solution of the system, determined the characteristics [32]. We write down the solution

, for ,, , (20)

where m(x,0) is the known initial conditions,

, for ,, , (21)

where φ(t) is the solution of Voltaire integral equation of the second kind [33]:

. (22)

**3.3.2. The second order asymptotic**

In order to find the second order asymptotic, we carry out some transformations in the basic equation (10).

In equation (16) by virtue of substitution u=εw, где ε=1/N, we return to the function *u(x)* and obtain

.

In equation (10) we replace

, (23)

and obtain

Denote this equation as (24).

From (8) and (23) follow, that H_{2}(u,t) is the characteristic functional for the quantity ξ(x,t)-Nm(x,t). The mathematical expectation of ξ(x,t)-Nm(x,t) equal to zero.

Denote ε^{2}=1/N and make replacements in the equation (24)

,, (25)

obtain a equality

(26)

This equality (26) is singularly perturbed equation [31]. We proceed similarly to the first order asymptotic. We consider a subclass of solutions H2(u,t) of equation (22), for which, by (25), the function F_{2}(w,t,ε) has the following properties.

There is a finite limit at ε→0 for F_{2}(w,t,ε)

, (27)

and its derivative with respect t and w(x)

, (28)

. (29)

We prove the following theorem.

**Theorem 2**. For ε→0 the solution of the problem (26) has the form

, (30)

where

, (31)

is cross-correlation function,

, , , , , (32)

are the initial and boundary conditions, and the summands of the cross-correlation function (31) are found from the system of partial differential equations

(33)

**Proof.** In equality** **(26) the exponential is expanded in series up to О(ε

^{3}). In view of (27-29), we make the passage to the limit for ε→0. We obtain that the functional F

_{2}(w,t) is the solution of a linear homogeneous equation

(34)

We write the solution of (34) in the form of a functional of the Gaussian distribution

, (35)

where R(y,z,t) – the cross-correlation function of stochastic density of the number of customers, serviced in the system at the moment t.

Substitute the solution (34) in equation (35), we obtain

(36)

Then we write the function R (y, z, t) in a symmetric form

, (37)

где δ(·) – δ-Dirac function, а σ^{2}(x,t), r(x,t), σ^{2}(t) are determined by the initial and boundary conditions

, , , , . (38)

We substitute the function R(y,z,t) (37) into the identity (36) (we assume that ), equate the coefficients for the same w^{2}(0), w^{2}(x) и w(0)w(x) (we remember (17), (38) and b(0,t)=0) and we get a system of three partial differential equations for functions σ^{2}(x,t), r(x,t) и σ^{2}(t), coincided with (33).

By virtue of (25) and also of (30), we can write the asymptotic equality for ε→0

.

Therefore, we can write

**Definition***.*** **A function

* *

is called the second order asymptotic of the characteristic functional F(u,t).

We write down the solution of the system (33)

, for , ,

, for , ,

, for , , , (39)

, for ,, .

The solutions of the systems (17) and (33) completely determine the parameters of the Gaussian distribution ξ(x,t). In this way the research of the autonomous non-Markov queuing system with an unlimited number of devices is completed.

**4. The application of the model and methods to the study of human population growth **

We consider a group of all women of the United State of America. The mathematical model as the autonomous non-Markov queuing system with an unlimited number of servers is defined by two functions S(x) and b(x,t). In demography these are a female survival function and a fertility function [6]. We define the survival function *S*(*x*) as a Gompertz-Meikem model [34]. And we write the function b(x,t) in the form

, (40)

where η(t) is the total fertility rate at time *t* and ψ(x,t) is the probability density of distribution for a reproductive age of the woman. We found the distribution of the probabilities of a reproductive age is invariant during the time, at least from 1980 to 2005, and defined the function ψ(x)=ψ(x,t) as the density of the bi-parametric γ-distribution

where , α>0.

Estimates of parameters α and β are found by the χ^{2 }– method.

Now we can make a forecast of the number of a female population in the United State of America, based on the assumptions about the dynamics of the total fertility rate and using the results (20-21). The initial conditions are the values of the number of women in 2010 [35, 36].

We assume η(*t*) will change from 1.93 to 2.15 in the period from 2010 to 2110. Applying formulas (20) and (21), we obtain the following values of the average number of women in the long run (Figure 3).

**Figure 3.** Scenario of the demographic situation in the period from 2010 to 2110, millions

Note that we can find the variance in the number of age groups, as well as the covariance between the number of the age group of zero age and the age group of age x for the period from 2010 to 2110, using the results of the second order asymptotic (39).

Therefore the proposed autonomous queuing system with an unlimited number of devices and the developed method of its investigation are an effective tool for estimating a long-term population size and for analyzing the demographic situation that has developed [36, 37].

**References**

- Kendall DG. Stochastic processes and population growth. J. Roy. Statist. Soc. 1949; ser. B, vol. 11: 230-264.
- Shripad D. Tuljapurkar, Hal Caswell. Structured-Population Models in Marine, Terrestrial, and Freshwater Systems [Internet]. Chapman & Hall; 1997 [cited 2017 Oct 17]. Available from https://link.springer.com/book/10.1007%2F978-1-4615-5973-3.
- Donald TR. Demographic methods and concepts. Oxford: Oxford University Press; 2006. 523 p.
- Hinde A. Demographic methods. Arnold; 1998. 305 p.
- Newell C. Methods and models in demography. Belhaven; 1988. 217 p.
- Vishnevskii AG. Demographic encyclopedic dictionary. Moscow: Sovetskaya entsiklopediya Publ.; 1985. 608 p.
- Coale А, Trussell J. The development and use of demographic models. Population Studies. 1996; Vol. 50, N 3: 469-484.
- Lee R. Probabilistic approaches to population forecasting. Population and Development Review. 1998; Vol. 24, Supplement: Frontiers of Population Forecasting: 156-190.
- Ahlburg D, Vaupel J. Alternative projections of the U.S. population. Demography. 1990; Vol. 27, N 4: 639-652.
- Sykes Z, Kim J. Dynamics of some special populations with NRR = 1. Demography. 1978; Vol. 15, N 4: 559-569.
- Sykes Z. Some stochastic versions of the matrix model for population dynamics. Journal of the American Statistical Association. 1969; Vol. 64, N 325: 111-130.
- Alho M, Spencer B. Uncertain population forecasting. Journal of the American Statistical Association. 1985; Vol. 80, N 390: 306-314.
- Pollard J. Continuous–time and discrete–time models of population growth. Journal of the Royal Statistical Society. 1969; Series A, Vol. 132, N 1: 93-97.
- Cohen J. Ergodicity of age structure in populations with Markovian vital rates, III: Finite state moments. Advances in Applied Probability. 1977; Vol. 9, N 3: 462–475.
- Keyfitz N. Introduction to the mathematics of population with revisions. Wesley Pub. Co.; 1977. 490 p.
- Caswell H. Matrix population models: construction, analysis, and interpretation. Sunderland, Massachusetts, Sinauer Associates; 1989. 722 p.
- Goodman L. Stochastic models for the population growth of the sexes growth. Biometrika. 1968; Vol. 55, N 3: 230-264.
- Booth H. Demographic forecasting: 1980 to 2005 in review. International Journal of Forecasting. 2006; N 22: 547-581.
- Reed L, Pearl L. On the summation of logistic curves. Journal of the Royal Statistical Society. 1927; Vol. 90, N 4: 729-746.
- Kapitza SP. The general theory of the growth of humankind. How many people have lived, live and are to live on Earth. Moscow: Nauka; 1999. 120 p.
- Kapitsa SP. Synergetics and forecasts of the future. Moscow: Editorial URSS; 2003. 283 p.
- Staroverov OV. Elements of mathematical demography. Moscow: Nauka; 1997. 158 p.
- Staroverov OV. Models of population movement. Moscow: Nauka; 1979. 230 p.
- McDonald JF, South DW. A comparison of two methods to project regional and state population for the U.S. The annals of regional science. 1973; Vol. 19, N 3: 40-53.
- Smith SK. State and local population projection. Springer Netherlands; 2001. 426 p.
- Whelpton PK. An empirical method of calculating future population. Journal of the American Statistical Association. 1936; Vol. 31, N 195: 457-473.
- Gnedenko BV, Kovalenko IN. Introduction to Queueing Theory. Moscow: Editorial URSS; 2005. 397 p.
- Nazarov AA, Terpugov AF. Queueing theory: educational material. Tomsk: NTL; 2004. 228 p.
- Kovalenko N, Filippova AA. Probability theory and mathematical statistics. Moscow: Vyssh. Shkola; 1973. 206 p.
- Nazarov AA Method of asymptotic analysis in queueing theory. Tomsk: NTL; 2004. 112 p.
- Lomov SA. Introduction to the general theory of singular perturbations. Moscow: Nauka; 1981. 400 p.
- Elsgolts L. Differential equations and the calculus of variations. Moscow: Nauka; 1969. 424 p.
- Krasnov ML. Integral equations. Moscow: Nauka; 1969. 216 p.
- Falin GI. Introduction to Actuarial Mathematics. Moscow University Publishing; 1994. 86 p.
- U.S. Population 2010 by age and gender [Internet]. U.S. Census Bureau; [cited 2017 Oct 17]. Available from http://proximityone.com/stateage.htm.
- Nosova MG. Autonomous non-Markov system of queuing system and its application in demography: Dis. ... cand. fiz.-mat. sciences: 05.13.18 / Nosova Mariia Gennadyevna. Tomsk; 2010. 204 p.
- Nazarov AA, Nosova MG. Mathematical model of the process of changing the demographic situation and its investigation. Reports of Tomsk State University of Control Systems and Radioelectronics. 2009; T. 2 (20): 100-105.