Passage 24 A #3

How does a population grow in size? - A mathematical approach

Geometric growth

Let (Subscript[N, t]) represent the population size at time t. If the population increases by 1% at each time, the size at time t+1, (Subscript[N, t+1]), is given by

	(Subscript[N, t+1])= (Subscript[N, t])+0.01 (Subscript[N, t]) = 1.01 (Subscript[N, t])			(1)

This is a difference equation that relates (Subscript[N, t]) to (Subscript[N, t+1]). We know that the solution is given by

	(Subscript[N, t])= (Subscript[N, 0])(1.01^t)						(2)
	
using the initial population size at t = 0, (Subscript[N, 0]). This gives a “geometric” sequence with the common ratio 1.01. 

In general, if the population increases by a% at each time the size (Subscript[N, t+1]) follows

	(Subscript[N, t+1])= (Subscript[N, t])+a/100 (Subscript[N, t]) = (1+a/100) (Subscript[N, t])	(3)

and the solution is

	(Subscript[N, t])= (Subscript[N, 0])(r^t)							(4)

with the common ration r = 1 + a/100. We know that if a > 0 (and r > 1), (Subscript[N, t]) exponentially increases. If a < 0 (and 0 < r < 1), it exponentially decreases to zero.

Assume a population where each individual gives birth to a certain number of offspring and survives a unit time (one year) with a certain survival probability. We ignore sex (male and female). Just assume a population of such as bacteria that proliferate by dividing themselves into several new individuals per unit time.

Then the population size at time t+1, (Subscript[N, t+1]), is given by

	(Subscript[N, t+1])=b (Subscript[N, t])+s (Subscript[N, t])					(5)

where the first term in the right hand side represents the number of new offspring recruited into the population in the next time step and the second term represents the number of individuals who survive one unit time. 

b is the number of offspring that survives to the next year per individual and s is the survival probability per unit time. This difference equation can be written as

	(Subscript[N, t+1])=r (Subscript[N, t])							(6)

as r = b + s.

The above equation (5) can be rearranged as

	ΔN = (Subscript[N, t+1])-(Subscript[N, t])= b (Subscript[N, t]) + s (Subscript[N, t]) - (Subscript[N, t]) = b (Subscript[N, t])-(1-s)(Subscript[N, t])		(7)

Thus, the net change in the population size N between successive times, ΔN, is given by the increase by birth minus the decrease by death. Note that d = 1 - s represents the probability of death.

Usually, the birth rate b is measured as “crude birth rate” as the number of births per 1000 people and the death rate d as “crude death rate” as the number of deaths per 1000 people.  

Data of crude birth/death rate are available online, e.g., https://data.worldbank.org/indicator

Descriptions of Wikipedia give more details of birth and death (mortality) rate. 
https://en.wikipedia.org/wiki/Birth_rate  
https://en.wikipedia.org/wiki/Mortality_rate

ΔN being positive means that the population increases, and it being negative means it decreases. In other words, b - d > 0 means increase and b - d < 0 means decrease; the same condition as r = b + s  = b + 1 - d > 1 or r < 1.

Now let’s see how a geometric sequence looks like.

Example 0 : Exponential model

Let' s set r to be 1.05, i.e., the population increases by 5% every year (Annual growth rate is 5%).

fn[x_]:=r x/.r->1.05

fn[100]

data=NestList[fn,1,50]

ListPlot[data,Joined->True,PlotMarkers->Automatic,AxesLabel->{"t","N_t"}]

ListLogPlot[data,Joined->True,PlotMarkers->Automatic,AxesLabel->{"t","N_t in log scale"}]

IMPORTANT: Geometric population growth is translated to a straight line when plotted in the logarithmic scale against time t, or a linear function of t.

This can be shown as follows. Take logarithm of both the sides of eq. (3)

	log (Subscript[N, t])= log (Subscript[N, 0])(r^t) = log (Subscript[N, 0])+log (r^t)= log (Subscript[N, 0])+ t log r

Therefore, log (Subscript[N, t]) is a linear function of time t, i.e., a straight line with the slope log r when plotted against time t.

dataLog=Log[data]

ListPlot[dataLog,PlotMarkers->Automatic,Joined->True,AxesLabel->{"t","Log N_t"}]

If we know the slop of the straight line in the logarithmic scale, we can obtain the common ration r from the slope.

Fit[dataLog,{1,t},t]

Exp[0.0487902]

Log[1.05]

Example 1: Collared doves

The doves brought in England as a pet increased from 4 to 18,855 in 10 years.

Log of the population size almost linearly increases with the slope 0.978. This means that the doves increased nearly exponentially with the common ratio e^0.978 = 2.66 every year.

Example 2 : A pheasant population in an island

Lack (1954) observed the number of pheasants in an island as follows.

data={{1937,1938,1939,1940,1941,1942},{8,30,81,282,705,1325}}//Transpose

ListPlot[data,PlotMarkers->Automatic,AxesLabel->{"Year","Pop size"}]

ListLogPlot[data,PlotMarkers->Automatic]

In logarithmic scale, it almost linearly increases with time t. Let’s examine the slope of the line.

dataPopSize={8,30,81,282,705,1325}

dataPopSizeLog=Log[dataPopSize]

gLog=ListPlot[dataPopSizeLog,PlotMarkers->Automatic]

We can obtain the slope by "regression analysis", a method of statistics.

reg=Fit[dataPopSizeLog,{1,t},t]

gReg=Plot[reg,{t,1,6}]

Show[gLog,gReg]

Coefficient[reg,t]

Exp[Coefficient[reg,t]]

The pheasant in the island increased roughly exponentially with the common ratio 2.82.

Example 3 : Survival of robins

Lack (1965) observed how robins survived for several years and the result is shown as follows.

dataNumAlive=Transpose[{{0,1,2,3,4},{129,49,20,8,2}}]

ListPlot[dataNumAlive,PlotMarkers->Automatic,AxesLabel->{"Time (Years)","Num indivs"}]

dataPercentageAlive=Transpose[{{0,1,2,3,4},{129,49,20,8,2}/129}]

ListPlot[dataPercentageAlive,PlotMarkers->Automatic,AxesLabel->{"Time (Years)","Percentage Alive"}]

dataPercentageAliveLog=Transpose[{{0,1,2,3,4},Log[{129,49,20,8,2}/129]}]

g1=ListPlot[dataPercentageAliveLog,PlotMarkers->Automatic,AxesLabel->{"Time (Years)","Percentage Alive"}]

fitted=Fit[dataPercentageAliveLog,{1,t},t]

g2=Plot[fitted,{t,0,4}]

Show[g1,g2]

Exp[-1.01457]

Annual survival probability of robins is 0.36.

Doubling time

For an exponential population growth, we can define the time necessary for the population size to be doubled. This is called "doubling time". Doubling time T can be solved as follows.

	(Subscript[N, t+T])= (Subscript[N, 0])(r^(t+T))=2(Subscript[N, t])=2(Subscript[N, 0])(r^t)
	
From this, the doubling time T = log 2 / log r

When r = 1.05 as above, the double time T = log 2/log 1.05

Log[2]/Log[1.001]

It takes about 14 years for the exponential growht with common ration r = 1.05 per year (annual 5% increase).

Half life

For an exponential population decline we can define the time necessary for the population size to be halved. This is called “half life". Half life T can be solved as follows.

	(Subscript[N, t+T])= (Subscript[N, 0])(r^(t+T))=(Subscript[N, t])/2=(Subscript[N, 0])(r^t)/2
	
From this, half life T = -log 2 / log r

When r = 0.95, the half life T = -log 2/log 0.95

-Log[2]/Log[0.9]

It takes about 13.5 years for the exponential decline with common ration r = 0.95 per year (5% decline annually).

Questions

Annual growth rate of the world population in 2015 is estimated as 1.07 %. If this rate is kept constant, how long will it take for the world population to be doubled, e.g., from 7.3 billion to 14.6 billion?

CountryData["World","Population"]

Log[2]/Log[1.0107]

Annual growth rate of Bangladesh population in 2015 is estimated as 1.61 %. If this rate is kept constant, how long will it take for Bangladesh population to be doubled, e.g., from 156 million to 312 million?

CountryData["Bangladesh","Population"]

Log[2]/Log[1.0161]

Annual growth rate of Japanese population in 2015 is estimated as -0.17%. If this rate is kept constant, how long will it take for Japanese population to be halved e.g., from 126 million to 63 million?

CountryData["Japan","Population"]

-Log[2]/Log[0.9983]

Check your country if it has been increasing exponentially in population size

dataUSA=CountryData["UK",{"Population",All}]

DateListPlot[dataUSA,PlotMarkers->Automatic]

DateListLogPlot[dataUSA,PlotMarkers->Automatic,PlotRange->All]