stats report

Population Dynamics

Next, we will consider mathematical modeling of population dynamics. These models are heavily used in the science of ecology and help us understand how populations interact with their environment and change over time.

Take N(t) as the size of a bacterial population at time t. Suppose that the bacteria divide and produce new bacterial cells at rate r per unit time. Then the population size at time t + ∆t is

N(t + ∆t) = N(t) + rN(t)∆t

Or

N(t + ∆t) −N(t) ∆t

= rN(r)

We take the limit as ∆t → 0 and get

dN

dt = rN(t)

Next, we want to solve this model. Multiply both sides by dt and divide by N:

dN

N = rdt

Integrate both sides

∫ t 0

dN

N = ∫ t

0 rdt

or ln N(t)|t0 = rt

ln N(t) − ln N(0) = rt

ln N(t) = ln N(0) + rt

N(t) = N(0)ert

This describes the growth of a population over time. N(0) is the population size at the initial time 0, r is the population growth rate, and N(t) is the population at time t. Consider a population that starts at size 2 and has a growth rate of 0.69 per year (the amount needed to double once per year):

1

N0=2 r<-c(1) t=seq(from=0,to=100,by=0.1) par(mfrow=c(1,1)) matplot(t,N0*exp(t(outer(r,t))),type="l",xlab="time (years)",ylab="population size N(t)",col=1:3,lty=1:3, main=paste("exponential growth, r=",r))

0 20 40 60 80 100

0 e

+ 0

0 2

e +

4 3

4 e

+ 4

3

exponential growth, r= 1

time (years)

p o

p u

la tio

n s

iz e

N (t

)

#legend("topleft",c(paste("r=",r[1]),paste("r=",r[2]),paste("r=",r[3])),lty=1:3,col=1:3)

In 100 years, the population grows to size 1040. If this were a population of rabbits (who can easily reproduce much faster than this), the mass of 1040 rabbits would be more than a billion times the mass of the sun.

Of course, populations cannot really grow infinitely. Instead, populations are regulated by resource limitations.

Discussion Questions

1. What resources put limits on human population growth? How about trees?

2. Make a plot of how you would expect the population growth rate to vary with population size.

2

Logistic Growth Model

A more realistic model of population growth must include resource limitations. One commonly used model is the logistic growth model:

dN

dt = rN(1 −

N

K )

where

N is population size r is the maximum population growth rate K is the carrying capacity of the population. This is the maximum population that can be sustained given the available resources.

The population growth rate in this model is

growth rate = r(1 − N

K )

When the population size N is small compared to the carrying capacity K, then we have N K ≈ 0 and the

growth rate is approximately r:

dN

dt ≈ rN

In this case, the population grows approximately exponentially. As the population size increases, then the population growth rate r(1 − N

K ) decreases.

r=1 K=1000 N=1:2000 growthrate=r*(1-N/K)

plot(N,growthrate,type="l",xlab="population size N", ylab="population growth rate", main="population growth rate under logistic growth model")

3

0 500 1000 1500 2000

− 1

.0 −

0 .5

0 .0

0 .5

1 .0

population growth rate under logistic growth model

population size N

p o

p u

la tio

n g

ro w

th r

a te

We see that the population growth rate decreases linearly as the size of the population increases.

Discussion Question Make a plot of how you think that population will vary with time under the logistic growth model starting at a low population size.

4

Solving the Logistic Growth Model

Unlike most population models, the logistic growth model actually has a solution:

N = K

1 + K−N0 N0

e−rt

The plot below r=0.69 K=1000 t=seq(from=0,to=20,by=0.1)

Nt=K/(1+(K-N0)/N0*exp(-r*t))

plot(t,Nt,type="l",xlab="time (years)", ylab="population size N", main="logistic growth model with K=1000")

0 5 10 15 20

0 2

0 0

4 0

0 6

0 0

8 0

0 1

0 0

0

logistic growth model with K=1000

time (years)

p o

p u

la tio

n s

iz e

N

We see that the population grows smoothly to the carrying capacity. At this point, the population growth rate is zero and so the population remains at the carrying capacity. At this point, birth and death are balanced and the population is stably living within the available resources.

Let’s also solve this using the ODE solver: library(deSolve)

## Warning: package 'deSolve' was built under R version 4.0.5

5

parameters<-c(r=0.69,K=1000) #specify model parameters state<-c(N=1) #initial state of the population

LogisticGrowthmodel<-function(t,state,parameters) {

with(as.list(c(state,parameters)),{

dN=r*N*(1-N/K) return(list(c(dN)))

}) }

times<-seq(from=0,to=30,by=0.01)

logGrowthout<-ode(y=state,times=times,func=LogisticGrowthmodel,parms=parameters)

par(oma=c(0,0,3,0)) plot(logGrowthout,xlab="time",ylab="population size N", main="logistc growth model from ODE solver")

0 5 10 15 20 25 30

0 2

0 0

6 0

0 1

0 0

0

logistc growth model from ODE solver

time

p o

p u

la tio

n s

iz e

N

We see that the dynamics look identical to the explicit solution, as we expect.

6

Human Carrying Capacity

The carrying capacity is the maximum population level that can be sustained by available resources.

If you google “human carrying capacity”, you can find many interesting resources. A good book on this topic is How Many People Can the Earth Support by Joel Cohen. Here is an excerpt from Wikipedia:

Several estimates of the carrying capacity have been made with a wide range of population numbers. A 2001 UN report said that two-thirds of the estimates fall in the range of 4 billion to 16 billion with unspecified standard errors, with a median of about 10 billion.[5] More recent estimates are much lower, particularly if non-renewable resource depletion and increased consumption are considered.[6][7] Changes in habitat quality or human behavior at any time might increase or reduce carrying capacity. In the view of Paul and Anne Ehrlich, “for earth as a whole (including those parts of it we call Australia and the United States), human beings are far above carrying capacity today.”

The application of the concept of carrying capacity for the human population has been crit- icized for not successfully capturing the multi-layered processes between humans and the environment, which have a nature of fluidity and non-equilibrium, and for sometimes being employed in a blame-the-victim framework. Supporters of the concept argue that the idea of a limited carrying capacity is just as valid applied to humans as when applied to any other species. Animal population size, living standards, and resource depletion vary, but the concept of carrying capacity still applies. The number of people is not the only factor in the carrying capacity of Earth. Waste and over-consumption, especially by wealthy and near-wealthy peo- ple and nations, are also putting significant strain on the environment together with human overpopulation. Population and consumption together appear to be at the core of many hu- man problems.Some of these issues have been studied by computer simulation models such as World. When scientists talk of global change today, they are usually referring to human-caused changes in the environment of sufficient magnitude eventually to reduce the carrying capacity of much of Earth (as opposed to local or regional areas) to support organisms, especially Homo sapiens.

The carrying capacity is the population size at which the population is in balance with its resources and thus can be sustained indefinitely. For humans, this is not a fixed value. It depends on technology, life styles, etc. Consider the carrying capacity for the contemporary US versus England in the Middle Ages. We in the US consume vastly more resources than someone living in the middle ages did. At the same time, our technology and economic system is far more efficient and advanced. We are able to grow vastly more food on the same amount of land and we are extremely good at extracting energy from our enivronment. Given this, our currenty carrying capacity appears to be much higher than what it would have been 1000 years ago.

The human carrying capacity of around 10 billion cited in the Wikipedia entry above is primarly based on food production. However, remember that carrying capacity means the population number that can be sustained indefinitely. Our current prosper is based heavily on energy production. Enery production is based heavily on fossil fuels. Fossil fuels are a non-renewable resource that is being exhausted quickly. The fossil fuel supply may last 50 years or it may last 200 years, but will definitely run out. Our current farming practices are energy dependent. Thus, much depends on whether we can transition to renewable energy recourses.

7

Stabiltiy Analysis of the Logistic Growth Model

We saw when used the ODE solver with the logistic growth that the population went smoothly to the equilibrium with is available resources. Will this always happen?

Consider the model again:

dN

dt = rN(1 −

N

K )

We see that this has a steady states atN = 0 and N = K:

dN

dt = 0 ⇒ rN̂(1 −

N̂

K ) = 0 ⇒ N̂ = 0 or N̂ = K

If the population moves slightly above N̂ = 0, then

dN

dt > 0

Thus, the population increases away from 0. This is an unstable steady state.

If the population is slightly below the carrying capacity, then (1 − N K

) > 0 and dN dt

> 0.

If the population is slightly below the carrying capacity, then (1 − N K

) < 0 and dN dt

< 0. That is, if we move slightly away from the steady state in either direction, we move back to the steady state. Thus, it is a stable equilibrium.

The condition for stability is that the system return to the steady state point when we move away. That is, dN dt

> 0 for N < N̂ dN dt

< 0 for N > N̂

In other words, the condition for stability is that the slope of the right hand side of dN dt

is negative at N̂:

if dN dt

= f(N), then the steady state N̂ is stable if df dN

< 0.

In the case of the logistic growth model, we have

f(N) = rN ( 1 −

N

K

) df

dN = r ( 1 −

N

K

) −rN

( 1 K

) = r ( 1 −

2N K

) At N̂ = 0, we have

df

dN = r ( 1 −

2 ∗ 0 K

) = r > 0

This is always positive and so the steady state at N̂ is unstable. This is simple to understand. When N = 0, then the population is at a steady state because there is no one in the population and thus it can’t grow. If move away from the steady state, this means that there are now some individuals in the population. Because the population is low, then the reproductive rate is high and the population grows exponentially away from N = 0.

At N̂ = K, we have df

dN = r ( 1 −

2 ∗K K

) = −r < 0

8

This is always negative and thus the steady state is always stable.

Discussion Question

1. Under the logistic growth model, the population always goes smoothly to its carrying capacity. Do you think that it this is realistic? That is, will populations always grow smoothly to equilibrium with their resources?

2. Do you think that the model is a realistic depiction of a growing population? List the assumptions of the model. What shortcomings does it have?

9

Discrete Time Models

In order to better understand what is going on, we will introduce discrete time model. Humans have overlapping generations. That is, there are people many different ages alive and reproducing at the same time. Thus, population growth is a continuous process. In contrast, many species have approximately non-overlapping generations. For example, in many insect species, most of the population emerges over a short time frame in the spring time, reproduces,and then dies. This produces a new generation, which are all born around the same, reproduce all at about the same time, and then die. This cycle continues until winter. For such species, it is reasonable to view time as progressing in discrete increments t = 1, 2, 3, ....

The discrete time logistic growth model is

N(t + 1) = rN(t) ( 1 −

N(t) K

) That is, each individual in generation t produces r

( 1 − N(t)

K

) offspring.

Discrete time models are easier to handle computationally than continous time ones: endtime=100 #length of time to simulate N<-vector(length=endtime) N[1]=20 r=1.2 K=1000

for(gen in 2:endtime) #loop over generations {N[gen]=r*N[gen-1]*(1-N[gen-1]/K) }

plot(1:endtime,N,type="l",xlab="generations", ylab="population size",main="discrete logistic growth")

10

0 20 40 60 80 100

5 0

1 0

0 1

5 0

discrete logistic growth

generations

p o

p u

la tio

n s

iz e

Long-term Behavior of the Discrete Logistic Model

In the case above, with the population growth rate set to 1.2 (meaning each individual produces an average of 1.2 offspring in the next generation), the population goes smoothly to an equilibrium.

The equilibrium occurs when N(t+ 1) = N(t) - that is, the population size doesn’t change between generations t and t + 1. This gives

N(t) = rN(t) ( 1 −

N(t) K

) This is satisfied by

N̂ = 0

and

1 = r ( 1 −

N̂

K

) or

N̂ = K(1 − 1 r

)

For the above example, we have r = 1.2 and K = 1000. This yields

11

N̂ = 1000(1 − 1

1.1 ) = 166.67

If the population growth rate r < 1, then we have

N̂ = K(1 − 1 r

) < 0

This steady state does not exist if r < 1. In this case the population dies out because reproduction is below replacement rate. endtime=100 #length of time to simulate N<-vector(length=endtime) N[1]=20 r=0.5 K=1000

for(gen in 2:endtime) #loop over generations {N[gen]=r*N[gen-1]*(1-N[gen-1]/K) }

plot(1:endtime,N,type="l",xlab="generations", ylab="population size",main="discrete logistic growth")

0 20 40 60 80 100

0 5

1 0

1 5

2 0

discrete logistic growth

generations

p o

p u

la tio

n s

iz e

If the population growth rate is between 2 and 3, there are dampled oscillations to the equilibrium:

12

endtime=100 #length of time to simulate N<-vector(length=endtime) N[1]=20 r=2.9 K=1000

for(gen in 2:endtime) #loop over generations {N[gen]=r*N[gen-1]*(1-N[gen-1]/K) }

plot(1:endtime,N,type="l",xlab="generations", ylab="population size",main=paste("discrete logistic growth, r=",r))

0 20 40 60 80 100

0 1

0 0

3 0

0 5

0 0

7 0

0

discrete logistic growth, r= 2.9

generations

p o

p u

la tio

n s

iz e

If the population growth rate is above 3, we get cycles: endtime=100 #length of time to simulate N<-vector(length=endtime) N[1]=20 r=3.2 K=1000

for(gen in 2:endtime) #loop over generations {N[gen]=r*N[gen-1]*(1-N[gen-1]/K) }

plot(1:endtime,N,type="l",xlab="generations", ylab="population size",main=paste("discrete logistic growth, r=",r))

13

0 20 40 60 80 100

0 2

0 0

4 0

0 6

0 0

8 0

0 discrete logistic growth, r= 3.2

generations

p o

p u

la tio

n s

iz e

endtime=100 #length of time to simulate N<-vector(length=endtime) N[1]=20 r=4 K=1000

for(gen in 2:endtime) #loop over generations {N[gen]=r*N[gen-1]*(1-N[gen-1]/K) }

plot(1:endtime,N,type="l",xlab="generations", ylab="population size",main=paste("discrete logistic growth, r=",r))

14

0 20 40 60 80 100

0 2

0 0

4 0

0 6

0 0

8 0

0 1

0 0

0 discrete logistic growth, r= 4

generations

p o

p u

la tio

n s

iz e

Time Delay and Overshooting the Resources

We see that as the growth rate gets higher, there are very large cycles. The population size goes far above its carrying capacity and then crashes to very low levels. This happens because the population grows past its available resources.

In the above run, the equilibrium is at

N̂ = 1000(1 − 1 4

) = 750

When the population reaches levels approaching 1000, then it is far past the sustainable level. Starvation (our some other consequence, depending on the limiting resource) ensues and the population crashes.

When the population growth rate is higher, then the population is able to grow further beyond the resource before consequences set in.

The reason why the discrete model displays this characteristic and the continuous one doesn’t is that in the continious model the population growth rate responds instantaneously to the current resource level. In the discrete model, there is a time delay between the population growth and the consequence of the population growth.

Consider generations 22 to 30 in the above run N[22:30]

## [1] 967.33704 126.38436 441.64542 986.37897 53.74198 203.41513 648.14966 ## [8] 912.20671 320.34250

15

At generation 23, the population size is 126. It grows to 441 in the next generation and is still below carrying capacity.Its growth rate is

r ( 1 −

N

K ) = 4(1 −

441 1000

) = 2.236

It grows to 986 in the next generation and is now above carrying capacity. Now the growth rate is

r ( 1 −

N(t) K

) = 4 ( 1 −

986 1000

) = 0.056

That is, the population drops to 54 (5.6% of its current level) because it overshot its resources.

Discussion Question

We have seen that delay in response to resources leads to unstability cycles in population models.

1. What would cycles like the last plot mean in the context of human populations?

2. Do you think that there is time delay in human response to resources?

16

Time Delay in continuous Models

Consider the situation in generation 24 in the example above. At this point, the population is well below its carrying capacity and resources are plentiful. For this reason, everyone in the population reproduces at a high rate. Unfortunately, this leaves the next generation well above carrying capacity and so there is mass starvation and population crashes.

This does not happen in the continuous logistic model because reproduction occurs continuosly and responds instantaenously to the current resources. It is possible to build a time delay in resource reponse into a continuous model. For example

dN

dt = rN(t)

( 1 −

N(t− τ) K

) Thus, the current population growth rate depends on what the population was at τ time units in the past. Typcially, current reproduction does not depend only on the current state. Rather, it depends on the whole lifetime experience of the organism. For an animal, current reproduction may depend on its overall state of health and its available energy reserves for making offspring (reproduction takes a lot of energy for the female).

In the case of humans, we are consuming resources (e.g. energy, fresh water, soil quality) that we are barely aware. Their depeletion will probably not affect our behavior until a crisis stage is reached. Thus, there is clearly a time delay in our response to resourcesm, although much more complex than the simple model above.

Time delay models are beyond the scope of this course. We can get the essential ideas from discrete models much more simply.

17

Chaos

Note that as the population growth rate increases, the dynamics get increasingly complicated:

0 10 20 30 40 50

0 3

0 0

7 0

0

discrete logistic growth, r= 3

generations

p o

p u

la tio

n s

iz e

0 10 20 30 40 50

0 4

0 0

discrete logistic growth, r= 3.3

generations

p o

p u

la tio

n s

iz e

0 10 20 30 40 50

0 4

0 0

discrete logistic growth, r= 3.6

generations

p o

p u

la tio

n s

iz e

0 10 20 30 40 50

0 4

0 0

1 0

0 0

discrete logistic growth, r= 3.9

generations

p o

p u

la tio

n s

iz e

At r = 3.3, we have regular cycles. At r = 3.6, we can discern periodic behavior, but it is complicated. There are four peaks in each cycle. At r = 3.9 there is not any discernable pattern. The behavior looks random. Let’s look more closely:

18

0 10 20 30 40 50

0 2

0 0

4 0

0 6

0 0

8 0

0 1

0 0

0 discrete logistic growth, r= 3.9

generations

p o

p u

la tio

n s

iz e

There is no periodic (repeating) dynamics. Instead, the dynamics look random. This behavior is known as chaotic dynamics. Even though the system is completely deterministic (i.e. no randomness in the data), it looks like random fluctuations. Consider the next plot, in which I run the system three times with slightly different starting population sizes:

19

0 10 20 30 40 50

0 2

0 0

4 0

0 6

0 0

8 0

0 1

0 0

0 discrete logistic growth, r= 3.9

generations

p o

p u

la tio

n s

iz e

N0=10 N0=11 N0=12

The only difference in the three lines is the starting population sizes of 10,11, and 12. We see that the population dynamics are in sync at first, but quickly diverge. If we look at any point beyond about 20 generations, the populations are at radically different sizes even though they started out very similar.

Contrast this to the case with r = 3.3:

20

0 10 20 30 40 50

0 2

0 0

4 0

0 6

0 0

8 0

0 discrete logistic growth, r= 3.3

generations

p o

p u

la tio

n s

iz e

N0=10 N0=11 N0=12

In this case the populations syncronize and after 20 generations that have identical dynamics.

The Butterfly effect

With r = 3.3, differences in the initial conditions dissapear with time.

With r = 3.9, differences in the initial conditions grow with time

This sensitivity to initial conditions is a hallmark of chaotic dynamics. With chaotic dynamics, tiny differences in initial conditions lead to large differences later in. This is often called the butterfly effect: “a butterfly flapping its wings in Chicago leads to a cyclone in the Indian Ocean”. This is an exagerration, but the basic idea is true: in chaotic systems, tiny effects can lead to big differences in the outcome. These effects are random for practical purposes and it is very hard to tell random dynamics from chaotic ones.

Chaotic dynamics are characterized by large, random-seeming fluctuations. The phenomenon was first observed in models of weather and was later found to exist in many types of dynamical systems. In population models, it tends too occur when population growth rates are large. Large population growth rates lead to large population fluctuations and large fluctuations lead to chaos.

Does Chaos Happen in Real Populations?

Chaos readily occurs in models, but does it occur in real populations? There has been a great deal of interest in this question. Consider Figures 1 and 2 below. These show scatter plots of estimates of the Lyapunov number for a selection of laboratory (Figure 1) and field (Figure 2) populations.

21

Figure 1: Figure from Ellner and Turchin, American Naturalist, 1995, vol 145(3)

The only thing that you need to know about the Lyapunov number is

γ < 1 ⇒ non-chaotic dynamics

γ > 1 ⇒ chaotic dyamics

We see from the figures that very few populations appear to exhibit chaotic dyamics. However, many populations are right at the boundary of chaos.

Evolution of Population Dynamics

Why should it be that many populations are at the boundary of chaos? The answer to this is unknown. However, we can speculate as to the reason.

Evolution works on reproduction. When a trait leads to higher survival and reproduction, then it is favored by evolution and tends to spread in a population. All else being equal, it is beneficial to an organism to produce more offspring. Of course “all else” is usually not equal and there are complex tradeoffs in determining the optimal strategy for an organism in terms of reproduction. Still, producing more offspring will often be a beneficial strategy and thus evolution will often favor it.

22

Figure 2: Figure from Ellner and Turchin, American Naturalist, 1995, vol 145(3)

23

This leads to a problem: Evolution works on individuals. An individual’s evolutionary interests are often best served by maximizing the number of offspring. However, high reproduction leads to unstable population dynamics, which is bad for the population.

Large fluctuations in population numbers can easily lead to population extinction. When the population number gets very low, some chance event (disease, bad weather, etc) can kill enough of the remaining members that the population goes extinct. Thus, natural selection on individuals may lead to extinction of the population.

This leads to multiple hypotheses about why we see many populations at the boundary of chaos, but few over:

1. Populations that cross the boundary into chaos tend to become extinct. Thus, the distribution of Lyapunov numbers that we observe is the result of a distribution of optimal individual reproductive numbers. However, the high end of the distribution gets truncated by population extinction.

2. Evolution favors being at the edge of chaos. That is, it might be most beneficial to have the highest reproductive rate that doesn’t lead to chaos. Thus, natural selection would favor reproductive rates right at the boundary. However, although this seems intuitively appealing, it has a major problem: Natural selection works on individuals, not groups. Suppose that we are in a situation where the average reproductive rate in the population is just below the threshold for chaos. An individual would still benefit by producing more offspring. If some individal has a mutation that leads to more offspring, then that mutation will tend to be represented more in future generations and it will thus spread and lead to higher reproductive rates. Even though the population as a whole will better off if each individual limits reproduction, individuals gain short-term benefit by “cheating” and reproducing more. Various researchers have suggested ways in which individual natural selection could lead to stable population dyanamics (my very first scientific publication was on this topic), but it remains an open question.

24

Dynamics of Multiple Populations: Predator-Prey Interactions.

Figure 3 shows a famous data set in population ecology. This shows the number of hare and lynx pelts traded by the Hudson Bay Company over a period of 100 years. We see a strong cyclic pattern in both populations, with a period of about 10 years. The lynx cycle lags behind the hare cycle by perhaps 1-2 years. The lynx is a major predator of the hare and the hare is a major food source for the lynx. The explanation for the observed pattern is that 1) in absence of many lynx, the hare population increases rapidly; 2) as the hare population increases, the lynx population also increases because there is plentiful food for the lynx; 3) As the lynx population begins to get large, they eat too many hare and the hare population collapses; 4) When the hare population collapses, the lynx population collapses soon after because there is not enough food; 5) Once the lynx population collapses, the hare population then starts to increase again. 6) Repeat. Ecologists are interested in understanding what drives these dynamics. There has been a large amount of work in developing mathematical models for population dynamics and then fitting these models to data.

Assumptions

Take X to represent the prey (e.g. hare) population number and Y to represent the predator (e.g. Lynx) population number.

We will assume 1. In absence of predators, the prey population follows the logistic growth model (and thus follows the assumptions of that model). 2. The prey is the only food source for the predator. Thus, in absence of the prey, the predator population decreases at rate r2. 3. The prey has only one predator. 4. The loss in prey biomass due to predation is given by α1XY1+βX . This will be explained below. 5. The gain in predator biomass due to predation is given by α2XY1+βX . 6. There is no randomness in the amount of reproduction, predation, etc. 7. The enivronment is constant. There are no changes in prey resource availability, mortality, etc. 8. There is no structure to either population. There is no variation due to age, sex, genetics, location, etc.

25

With these assumptions, we have

dX

dt = r1X

( 1 −

X

K

) −

α1XY

1 + βX

dY

dt = −r2Y +

α2XY

1 + βX

r1 ( 1 − X

K

) is the growth rate of the prey population in absence of the predator (logistic growth)

−α1XY1+βX is the decrease in the prey biomass due to predation

−r2Y is the is the (negative) growth rate of the predator population in absence of the prey (they die out with no food). α2XY 1+βX is the increase in predator biomass due their predation on prey.

Functional Response

The term α2XY1+βX is called the functional response. This determines the amount of prey eaten as a function of prey density.

The function that we use here is called a Holling Type II functional response. Holling was a scientist who did early work on the concept of functional response. Plot the function: maxX=100 #max number of pret eaten halfmax=50 #prey number at half max eaten

beta=1/halfmax alpha=beta*maxX

X=seq(from=0,to=1000,by=1) plot(X,alpha*X/(1+beta*X),type="l", ylab="prey eaten", xlab="prey number",main="Functional Reponse")

26

0 200 400 600 800 1000

0 2

0 4

0 6

0 8

0 Functional Reponse

prey number

p re

y e

a te

n

When the prey number is low, then the number eaten per predator increases linearly - that is if there is twice the prey, then twice the prey is eaten. As the number of prey keeps on increasing, then the number eaten per predator begins to flatten out. This is because the predators are getting enough to and so don’t eat more as the number of prey increases.

We assume that the functional reponse is the same between predators and prey except for the constants α1 and α2. The ratio α2α1 determines how efficiently prey biomass is converted into predator biomass).

Dynamics of the Model

Let’s run the model using the ODE solver: library(deSolve)

parameters<-c(r1=1,K=500,alpha1=0.01,beta=0.005,r2=1,alpha2=0.01) #specify model parameters state<-c(X=200,Y=10) #initial state of the population

PredatorPreymodel<-function(t,state,parameters) {

X=state[1] Y=state[2]

with(as.list(c(state,parameters)),{

dX=r1*X*(1-X/K)-alpha1*X*Y/(1+beta*X) dY=-r2*Y+alpha2*X*Y/(1+beta*X)

27

return(list(c(dX,dY)))

})

}

times<-seq(from=0,to=200,by=0.01)

PredPreyout<-ode(y=state,times=times,func=PredatorPreymodel,parms=parameters)

par(mfrow=c(1,1)) matplot(PredPreyout[,1],PredPreyout[,2:3],type="l",xlab="time", ylab="population numbers", main="Predator-Prey Model")

0 50 100 150 200

0 1

0 0

2 0

0 3

0 0

4 0

0

Predator−Prey Model

time

p o

p u

la tio

n n

u m

b e

rs

In this case, the predator and prey populations undergo damped oscillations to an equilibrium. If we increase the carrying capacity to K = 625, the we get cycles: library(deSolve)

parameters<-c(r1=1,K=625,alpha1=0.01,beta=0.005,r2=1,alpha2=0.01) #specify model parameters state<-c(X=200,Y=10) #initial state of the population

PredatorPreymodel<-function(t,state,parameters) {

28

X=state[1] Y=state[2]

with(as.list(c(state,parameters)),{

dX=r1*X*(1-X/K)-alpha1*X*Y/(1+beta*X) dY=-r2*Y+alpha2*X*Y/(1+beta*X)

return(list(c(dX,dY)))

})

}

times<-seq(from=0,to=200,by=0.01)

PredPreyout<-ode(y=state,times=times,func=PredatorPreymodel,parms=parameters)

par(mfrow=c(1,1)) matplot(PredPreyout[,1],PredPreyout[,2:3],type="l",xlab="time", ylab="population numbers", main="Predator-Prey Model")

0 50 100 150 200

0 1

0 0

2 0

0 3

0 0

4 0

0 5

0 0

Predator−Prey Model

time

p o

p u

la tio

n n

u m

b e

rs

This shows the phase plane analysis:

29

library(phaseR)

## Warning: package 'phaseR' was built under R version 4.0.5

## ------------------------------------------------------------------------------- ## phaseR: Phase plane analysis of one- and two-dimensional autonomous ODE systems ## ------------------------------------------------------------------------------- ## ## v.2.1: For an overview of the package's functionality enter: ?phaseR ## ## For news on the latest updates enter: news(package = "phaseR") lim=500

PredPreyflow<-flowField(PredatorPreymodel,xlim=c(0,lim), ylim=c(0,lim), parameters=parameters,add=FALSE,xlab="X",ylab="Y")

PredPreyNullclines<-nullclines(PredatorPreymodel,xlim=c(0,lim), ylim=c(0,lim), parameters=parameters,state.names = c("X","Y"))

PredPreyTrajectory<-trajectory(PredatorPreymodel,y0=matrix(c(25,200,150,150,450,200,150,50,400,400,200,25),nrow=6,byrow=T),parameters=parameters,tlim=c(0,20))

## Note: col has been reset as required

0 100 200 300 400 500

0 1

0 0

2 0

0 3

0 0

4 0

0 5

0 0

X

Y

X nullclines Y nullclines

We see that with these parameter values there is a single stable equilibrium. The system will osciallate around this equilibrium wherever it starts from (except if either population starts at zero).

Using models such as this one, we can explore how the dynamics change as we vary different aspects of the model. In particular, ecologists are interested in what characteristics will cause cycles versus more stable dynamics and how well the correspondance is between the models and real populations.

30

Dynamics of Large Scale Fishing

Now we will turn out attention to a more specialized example of predation: large scale fishing. Let’s start with the concept of a fishery. Here is a definition from https://www.sciencedaily.com/terms/fishery.htm

Fishery A fishery (plural: fisheries) is an organized effort by humans to catch fish or other aquatic species, an activity known as fishing. Generally, a fishery exists for the purpose of providing human food, although other aims are possible (such as sport or recreational fishing), or obtaining ornamental fish or fish products such as fish oil. Regardless of purpose, however, the term fishery generally refers to a fishing effort centered on either a particular ecoregion or a particular species or type of fish or aquatic animal, and usually fisheries are differentiated by both criteria. Examples would be the salmon fishery of Alaska, the cod fishery off the Lofoten islands or the tuna fishery of the Eastern Pacific. Most fisheries are marine, rather than freshwater; most marine fisheries are based near the coast.This is not only because harvesting from relatively shallow waters is easier than in the open ocean, but also because fish are much more abundant near the coastal shelf, due to coastal upwelling and the abundance of nutrients available there.

Many fisheries have collapsed due to overfishing.

Stock Peak Catch (year) Catch in 1981 Antartic Blue Whale 29,000 whales (1931) Nil Antartic fin whale 27,000 whales (1938) Nil Hokkaido herring 850,000 tons (1913) Nil Peruvian anchoveta 12.3 million tons(1970) 0.3 million tons SW African pilchard 1.4 million tons (1962) Nil Nort Sea herring 1.5 million tons (1962) negligible California sardine 640,000 tons (1936) Nil Georges Bank herring 374,000 tons (1968) Nil Japanse sardine 2.3 million tons (1939) 17,000 tons

When a fishery collapes, it is a disaster both environmentally (for the fish) and economically (for the fishermen and consumers of the fish)

Even when fisheries do not collapse, fishing is ofen on economically marginal occupation. In other words, most fishermen are poor.

Goal: Manage fisheries in way that is environmentally and economically sustainable, provides good work for fishermen, and provides a good source to consumers.

Mathematical models have played a large role in fisheries management.

###Modeling of Fisheries

31