264 KiB
jupyter nbconvert filename.ipynb --to latex%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from scipy import sparse
import pandas as pd
from IPython.display import display
eye = np.eye(4)
print(eye)
sparse_mtx = sparse.csr_matrix(eye)
print(sparse_mtx)
x = np.linspace(-10,10,100)
y = np.sin(x)
plt.plot(x,y,marker='x')
plt.show()
data = {'Name': ["John", "Anna", "Peter", "Linda"], 'Location': ["Nairobi", "Napoli", "London", "Buenos Aires"], 'Age':[51, 21, 34, 45]}
data_pandas = pd.DataFrame(data)
display(data_pandas)import numpy as np
import matplotlib.pyplot as plt
from scipy import sparse
import pandas as pd
from IPython.display import display
import mglearn
import sklearn
from sklearn.linear_model import LinearRegression
from sklearn.tree import DecisionTreeRegressor
x, y = mglearn.datasets.make_wave(n_samples=100)
line = np.linspace(-3,3,1000,endpoint=False).reshape(-1,1)
reg = DecisionTreeRegressor(min_samples_split=3).fit(x,y)
plt.plot(line, reg.predict(line), label="decision tree")
regline = LinearRegression().fit(x,y)
plt.plot(line, regline.predict(line), label= "Linear Rgression")
plt.show()import numpy as np
from matplotlib import pyplot as plt
# Load in data file
data = np.loadtxt('src/Hudson_Bay.csv', delimiter=',', skiprows=1)
# Make arrays containing x-axis and hares and lynx populations
year = data[:,0]
hares = data[:,1]
lynx = data[:,2]
plt.plot(year, hares ,'b-+', year, lynx, 'r-o')
plt.axis([1900,1920,0, 100.0])
plt.xlabel(r'Year')
plt.ylabel(r'Numbers of hares and lynx ')
plt.legend(('Hares','Lynx'), loc='upper right')
plt.title(r'Population of hares and lynx from 1900-1920 (x1000)}')
plt.savefig('Hudson_Bay_data.pdf')
plt.savefig('Hudson_Bay_data.png')
plt.show()Hares and lynx in Hudson bay from 1900 to 1920
Why now create a computer model for the hare and lynx populations?
We see from the plot that there are indeed fluctuations. We would like to create a mathematical model that explains these population fluctuations. Ecologists have predicted that in a simple predator-prey system that a rise in prey population is followed (with a lag) by a rise in the predator population. When the predator population is sufficiently high, then the prey population begins dropping. After the prey population falls, then the predator population falls, which allows the prey population to recover and complete one cycle of this interaction. Thus, we see that qualitatively oscillations occur. Can a mathematical model predict this? What causes cycles to slow or speed up? What affects the amplitude of the oscillation or do you expect to see the oscillations damp to a stable equilibrium? The models tend to ignore factors like climate and other complicating factors. How significant are these?
-
We see oscillations in the data
-
What causes cycles to slow or speed up?
-
What affects the amplitude of the oscillation or do you expect to see the oscillations damp to a stable equilibrium?
-
With a model we can better understand the data
-
More important: we can understand the ecology dynamics of predator-pray populations
The traditional (top-down) approach
The classical way (in all books) is to present the Lotka-Volterra equations:
import numpy as np
import matplotlib.pyplot as plt
def solver(m, H0, L0, dt, a, b, c, d, t0):
"""Solve the difference equations for H and L over m years
with time step dt (measured in years."""
num_intervals = int(m/float(dt))
t = np.linspace(t0, t0 + m, num_intervals+1)
H = np.zeros(t.size)
L = np.zeros(t.size)
print('Init:', H0, L0, dt)
H[0] = H0
L[0] = L0
for n in range(0, len(t)-1):
H[n+1] = H[n] + a*dt*H[n] - b*dt*H[n]*L[n]
L[n+1] = L[n] + d*dt*H[n]*L[n] - c*dt*L[n]
return H, L, t
# Load in data file
data = np.loadtxt('src/Hudson_Bay.csv', delimiter=',', skiprows=1)
# Make arrays containing x-axis and hares and lynx populations
t_e = data[:,0]
H_e = data[:,1]
L_e = data[:,2]
# Simulate using the model
H, L, t = solver(m=20, H0=34.91, L0=3.857, dt=0.1,
a=0.4807, b=0.02482, c=0.9272, d=0.02756,
t0=1900)
# Visualize simulations and data
plt.plot(t_e, H_e, 'b-+', t_e, L_e, 'r-o', t, H, 'm--', t, L, 'k--')
plt.xlabel('Year')
plt.ylabel('Numbers of hares and lynx')
plt.axis([1900, 1920, 0, 140])
plt.title(r'Population of hares and lynx 1900-1920 (x1000)')
plt.legend(('H_e', 'L_e', 'H', 'L'), loc='upper left')
plt.savefig('Hudson_Bay_sim.pdf')
plt.savefig('Hudson_Bay_sim.png')
plt.show()The plot
If we perform a least-square fitting, we can find optimal values for the parameters a, b, d, c. The optimal parameters are a=0.4807, b=0.02482, d=0.9272 and c=0.02756. These parameters result in a slightly modified initial conditions, namely H(0) = 34.91 and L(0)=3.857. With these parameters we are now ready to solve the equations and plot these data together with the experimental values.
Linear regression in Python
import numpy as np
import matplotlib.pyplot as plt
from IPython.display import display
import sklearn
from sklearn.linear_model import LinearRegression
from sklearn.tree import DecisionTreeRegressor
data = np.loadtxt('src/Hudson_Bay.csv', delimiter=',', skiprows=1)
x = data[:,0]
y = data[:,1]
line = np.linspace(1900,1920,1000,endpoint=False).reshape(-1,1)
reg = DecisionTreeRegressor(min_samples_split=3).fit(x.reshape(-1,1),y.reshape(-1,1))
plt.plot(line, reg.predict(line), label="decision tree")
regline = LinearRegression().fit(x.reshape(-1,1),y.reshape(-1,1))
plt.plot(line, regline.predict(line), label= "Linear Regression")
plt.plot(x, y, label= "Linear Regression")
plt.show()