3 - Competition continued#
Linear stability analysis#
Beyond this graphical method it would obviously be useful if we can formally characterise these equilibria in some way, and in particular assess when each will be an attracting end-point of the dynamics (the equilibrium is stable) and when not (the equilibrium is unstable). Much like we did with the logistic model, we can assess the local stability of steady states of a system by linearising around the equilibrium and finding the eigenvalues by studying the Jacobian matrix. Let us suppose that \(dR/dt=f(R,G)\) and \(dG/dt=g(R,G)\). Then the Jacobian matrix is made up of the partial derivatives of \(f\) and \(g\),
taken at the equilibrium values of \(R\) and \(G\). We can then use this Jacobian to find the eigenvalues at an equilibrium, which will tell us whether that point is stable or unstable. In particular, if,
\(\lambda_1<0\) and \(\lambda_2<0\), the equilibrium is stable;
\(\lambda_1>0\) and \(\lambda_2>0\), the equilibrium is unstable;
\(\lambda_1<0\) and \(\lambda_2>0\), the equilibrium is a saddle.
Since we are dealing with a 2x2 matrix, and since we only really care about the signs of the eigenvalues, we can usually determine stability just by finding the trace and determinant of our matrix:
\(tr(J) = \frac{\partial f}{\partial R} + \frac{\partial g}{\partial G}\),
\(det(J) = \left[\frac{\partial f}{\partial R}\right].\left[\frac{\partial g}{\partial G}\right]-\left[\frac{\partial f}{\partial G}\right].\left[\frac{\partial g}{\partial R}\right].\)
These relate to the eigenvalues through as \(tr(J)=\lambda_1+\lambda_2\) and \(det(J)=\lambda_1\lambda_2\). The signs of these two quantities combine to tell us about the stability of the equilibrium, which is best seen through the following diagram:
Figure: Stability of an equilibrium using the trace and determinant of the Jacobian.
Note
For a 2x2 Jacobian, if one of the off-diagonal entries is zero, the eigenvalues are simply the two entries in the main-diagonal.
Let’s write down the generic Jacobian for our system,
We can use the information from the Jacobians to derive the eigenvalues at each equilibria explicitly. In particular we can use the signs of the trace and determinant to establish the stability of each equilibrium.
The extinction equilibrium
At \((0,0)\), the Jacobian is just
which has \(T=a+c>0\), \(D=ac>0\) and \(T^2-4D = (a-c)^2 >0\), where \(T\) and \(D\) are the trace and determinant of the matrix. Thus \((0,0)\) is always an unstable node.
The two single-species equilibria
At \((a,0)\), we have
We could calculate the trace and determinant, or we can note that the 0 in the off-diagonal means we can just read off the eigenvalues as \(\lambda_1=-a\) and \(\lambda_2=c-ad\). Therefore, looking back at our phase portraits from last time, in cases (a) and (c), when \(c-ad<0\) (look at the horizontal axis), we have a stable node at this point. In cases (b) and (d), when \(c-ad>0\), we have a saddle point. These findings are in accord with the direction fields and our previous interpretation of the behaviour of the system.
Using a similar argument (or just the symmetry of the equations), we can show that \((0,c)\) is a saddle in cases (a) and (d), and a stable node in cases (b) and (c).
The coexistence equilibrium
For the final equilibrium, we could try and calculate the equilibrium densities explicitly, but it is no point making unnecessary work for ourselves. Instead, note that it lies at the intersection of the two off-axis nullclines. Substituting the nullcline equations into \(\mathbf{J}\), we get
which has \(T=-(R_*+G_*)<0\), \(D=(1-bd)R_*G_*\) and \(T^2-4D=(R_*-G_*)^2+4bdR_*G_*>0\). Everything hinges on the sign of \(D=(1-bd)R_*G_*\).
In case (c), we have \(a-bc>0\) and \(c-ad>0\). Putting these together we find we have \(a > bad\) and so \(1>bd\). Therefore \(D>0\) and C is a stable node.
In contrast, in case (d) we have \(a<bc\) and \(c<ad\), meaning \(a < bad\) and so \(1<bd\). Therefore \(D<0\) and C is a saddle point.
Have a go: coding
To use the code, click the rocket icon at the top of the page, and then click Colab. This will take you to a live version of this page where you can run the code.
Scroll down to the code cell below and either click in the cell and hit Shift-Enter, or click on the ‘Play’ icon in the top-left of the cell. You may see a warning about loading content from Github - click ‘Run anyway’.
This code will produce a time-course and phase portrait for the dynamics of two competing species, here assumed to be Red and Grey Squirrels. The dynamics are governed by the following ordinary differential equations:
After it has run, below the code cell you will see a set of sliders to choose the values of the 4 parameters, and boxes for you to enter two lots of initial conditions for \(R\) and \(G\). When you have chosen your values you will see two plots appear - a time course for Red and Grey squirrels, and a phase portrait (including nullclines (red/black) and 2 trajectories (blue)). You can then change the parameter values or initial conditions.
Note that the sliders/plot in this webpage are static and do not change - you need to load the Colab version for them to work.
Show code cell source
################
# HAVE A GO: CODING
################
# Import the necessary libraries
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
import ipywidgets as widgets
from ipywidgets import interactive
aw=widgets.FloatSlider(min=0.5, max=5, step=0.25, value=3, description='a')
bw=widgets.FloatSlider(min=0.5, max=5, step=0.25, value=0.5, description ='b')
cw=widgets.FloatSlider(min=0.5, max=5, step=0.25, value=3,description='c')
dw=widgets.FloatSlider(min=0.5, max=5, step=0.25, value=0.75,description='d')
R0_1w=widgets.BoundedFloatText(value=3,min=0,max=5,description='1st R(0)')
G0_1w=widgets.BoundedFloatText(value=3,min=0,max=5,description='1st G(0)')
R0_2w=widgets.BoundedFloatText(value=1,min=0,max=5,description='2nd R(0)')
G0_2w=widgets.BoundedFloatText(value=1,min=0,max=5,description='2nd G(0)')
u1 = widgets.HBox([R0_1w, G0_1w])
u2 = widgets.HBox([R0_2w, G0_2w])
#Lotka-Volterra dynamics
def competition(N,t,a,b,c,d):
R=N[0]
G=N[1]
dR = R*(a - R - b*G)
dG = G*(c - G - d*R)
return [dR,dG]
def runsquirrels(a,b,c,d,R0_1,G0_1,R0_2,G0_2):
N0=[R0_1,G0_1]
N1=[R0_2,G0_2]
tc = np.linspace(0, 10, 1000)
Nc = odeint(competition, N0, tc,args=(a,b,c,d))
Nd = odeint(competition, N1, tc,args=(a,b,c,d))
plt.rcParams['figure.figsize'] = [12, 4]#
plt.rcParams.update({'font.size': 16})
fig, (ax1, ax2) = plt.subplots(1, 2)
ax1.plot(tc, Nc[:,0], "r", label="Reds")
ax1.plot(tc, Nc[:,1], "k", label="Greys")
ax1.plot(tc, Nd[:,0], "r:")
ax1.plot(tc, Nd[:,1], "k:")
ax1.set(xlabel='Time', ylabel='Densities')
ax1.legend()
ax1.axis([0,10,0,5])
rr=np.linspace(0,10,5)
rnull=(a-rr)/b
gnull=c-d*rr
ax2.plot(Nc[:,0],Nc[:,1],'b')
ax2.plot(Nd[:,0],Nd[:,1],'b')
ax2.plot(rr,rnull,'r')
ax2.plot(rr,gnull,'k')
ax2.axis([0, 5, 0, 5])
ax2.set(xlabel='Red Squirrels', ylabel='Grey Squirrels')
plt.show()
return()
out = widgets.interactive_output(runsquirrels, {'a': aw, 'b': bw, 'c': cw, 'd': dw, 'R0_1': R0_1w, 'G0_1': G0_1w, 'R0_2': R0_2w, 'G0_2': G0_2w})
display(aw,bw,cw,dw)
display(u1,u2,out)
A bifurcation diagram#
Transcritical bifurcations in the competition model. Parameter values: \(b=0.5\), \(c=5\), \(d=0.75\).
We can also plot a bifurcation diagram of the system. This is slightly complicated since we now have two variables and multiple parameters, but we can just focus on how one parameter impacts one variable. Here we see three different transcritical bifurcations at different points as we vary the red squirrel carrying capacity, \(a\).
When \(a\) is low, the Grey-only equilibrium (grey line with \(R=0\)) is stable, the Red-only equilibrium (red line) is unstable and the coexistence equilibrium (blue) is negative (and unstable).
When \(a\) is intermediate, the coexistence equilibrium becomes positive and stable and the Grey-only equilibrium becomes unstable through a transcritical bifurcation. The Red-only equilibrium remains unstable.
When \(a\) is high, the coexistence equilibrium becomes unstable (and actually negative in \(G\), though this plot does not show that) and the Red-only equilibrium becomes stable through a transcritical bifurcation. The Grey-only equilibrium remains unstable.
(There is a third transcritical bifurcation at 0 where the Red-only and Grey-only equilibria collide, but this makes no real difference to the dynamics since we assume \(a>0\)).
Note what we have achieved here: we have mapped out all possible behaviours of the system without having to specify exactly what the values of the parameters are, and without having to calculate the location of the equilibrium point.
3 key points#
Linear stability analysis allows us to classify the behaviour of a system near an equilibrium.
Competing species can coexist with each other in the long term, provided the effects of interspecific competition are quite weak (in fact, if we assume equal carrying capacities, \(a=c\), we require both \(b<1\) and \(d<1\)).
Very often there will be competitive exclusion, where only one species can survive in the long term, whenever the effects of inter-specific competition is quite strong.