MATLAB Answers

Lotka-Volterra, can you estimate steady-state values from just exploring a plot?

5 views (last 30 days)
Kamila Bienkowska
Kamila Bienkowska on 11 Jul 2021
Edited: David Goodmanson on 16 Jul 2021
I am currently working on the Lotka-Volterra model for modelling hare-lynx populations over time.
I have a working function to model this with ode45, that is not a problem.
These are the differential equations that I am using:
% growth rate of the prey population (snowshoe hare)
% growth rate of the predator populations (canadian lynx)
But then, I also have historical data of number of lynx and hares from 1900 to 1920. I was asked to estimate a new k1 value (different from my initial model, where I know all, k1, k2 and k3), then estimate the steady-state values for X (hare) and Y (lynx) from the historical data. And then, based on the new k1 and the calculated steady-states, calculate new k2 and k3.
So I have the new k1, which is 1.3610. But I don't know how to estimate the steady-states here from just looking at a graph (this was suggested by one of the tutors on a forum), or without knowing the new k2 and k3 yet. If I knew the k2 and k3, that would be fairly easy, as I would only need to calculate X and Y values for which my differential equations are zero. I read that there is an option of plotting nullclines to find steady-states, but I don't know how to apply them with only knowing k1 and having the historical values of number of prey and predators.
The plots above are the historical data that I have imported and plotted. I would very much appreciate if someone could tell me how I can estimate the steady-states from these.

Answers (1)

David Goodmanson
David Goodmanson on 16 Jul 2021
Edited: David Goodmanson on 16 Jul 2021
Hi Kamila,
Here is one fairly simple way to find the parameters. The two equations are
dhare/dt = k1*hare - k2*hare*lynx
dlynx/dt = -k3*lynx + k2*hare*lynx
The steady state fixed point occurs when both time derivatives are zero in which case the fixed point is at [hare0,lynx0] and you can solve and get the standard result
hare0 = k3/k2 lynx0 = k1/k2
For Lotka-Volterra, the oscillations form a closed single-loop orbit. Fortunately for the Hudson Bay data someone chose to provide 21 years so there is an almost exactly closed two-loop curve. (The loops would be identical for perfect L-V behavior but of course this is real data).
Suppose 'hare' is a row vector of population by year. I added the first point to the end of the array with hare = [hare hare(1)] (similarly for lynx) so now if you plot hare vs lynx you will get a closed curve with two loops.
So far so good. Now suppose you time integrate both sides of the first L-V equation around an integral number of orbits, going from tstart to tfin. On the left hand side, integrating the derivative you get hare(tstart)-hare(tfin) = 0 since the start position and the end position are the same. On the right hand side you get
0 = k1*Ihare - k2*Iboth where Ihare = Int hare dt and Iboth = Int hare*lynx dt
Iboth/Ihare = k1/k2 = lynx0
and similarly for the second equation
Iboth/Ilynx = k3/k2 = hare0
For the integrations you can use t = 0:21 since you have one-year intervals and use the trapz function ,
Ihare = trapz(t,hare) (similary for lynx); Iboth = trapz(t,hare.*lynx)
I got lynx0 = 2.0883e+04 hare0 = 3.5292e+04
which, given the look of the data, ought to be rounded to something like 21000 and 35000. Once you have those, with a new k1 you can use the eqns above to solve for k2 and k3.
A long time ago there was a comic strip about Fat Freddy's Cat. One day Fat Freddy's Cat went out toward the woods and met a big kind of cat he had never seen before. He said, what kind of cat are you? And the cat said, "I'm a lynx". To which FFC said, "Oh, I always wondered where links sausage came from".

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!