<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Ecology Example
We use here the example of migratory birds in order to demonstrate the model fitting of this package
Step1:  Model
Here we use a model of migratory birds, were the rate of observation changes across time and the observed counts are negative-binomially
distributed. Associated parameters are
* peak rate of observations $m_0$
* aggregation parameter $k$
* where peak occurs in year $T_0$
Step2:  check rates
First let's check that the functional form of the rates are corrected. Let's choose $m_0=10.0$ and $T_0=0.2$
Step3:  It slightly underestimates heterogeneity, but is close for max rate
<ASSISTANT_TASK:>
Python Code:
%matplotlib inline
import ABCPRC as prc
import numpy as np
import scipy.stats as stats
import matplotlib.pyplot as plt
Explanation: Ecology Example
We use here the example of migratory birds in order to demonstrate the model fitting of this package
End of explanation
def ibm(*ps):
    m0,k = ps[0],ps[1]
    T0 = 0.5
    #measurements in regular increments throughout the year
    ms,ts = np.zeros(100),np.linspace(0,1,100)
    ms = (m0/2)*(np.sin(np.pi/2 + (ts-T0)*2*np.pi)+1)
    ps = k/(k+ms)
    outs = np.array([stats.nbinom(n=k,p=p).rvs() if p>0 else 0 for m,p in zip(ms,ps)])
    return outs
m0,k = 50.0,1.0
xs = ibm(m0,k)
plt.plot(xs,'ko')
Explanation: Model
Here we use a model of migratory birds, were the rate of observation changes across time and the observed counts are negative-binomially
distributed. Associated parameters are
* peak rate of observations $m_0$
* aggregation parameter $k$
* where peak occurs in year $T_0$
End of explanation
#measurements in regular increments throughout the year
ms,ts = np.zeros(100),np.linspace(0,1,100)
ms = (m0/2)*(np.sin(np.pi/2 + (ts-0.5)*2*np.pi)+1)
plt.plot(ts,ms)
plt.ylabel('rate')
plt.xlabel('time (yrs)')
m = prc.ABC()
priors = [stats.expon(scale=100.0).rvs,stats.expon(scale=0.5).rvs]
m.setup(modelFunc=ibm,xs=xs,priors=priors)
m.fit(sample_size=100)
m.run(1000)
res = m.trace()
plt.figure()
print('Initial Distribution')
m.trace(plot=True,tol=0)
plt.figure()
print('Middle Tolerance')
m.trace(plot=True,tol=5)
plt.figure()
print('Final Distribution')
m.trace(plot=True,tol=-1)
ps = np.round(m.paramMAP(),decimals=2)
print('MAP for max rate is : {}, MAP for heterogeneity is {}'.format(*ps))
res = m.fitSummary()
Explanation: check rates
First let's check that the functional form of the rates are corrected. Let's choose $m_0=10.0$ and $T_0=0.2$
End of explanation
m.save('ecology_model_example')
Explanation: It slightly underestimates heterogeneity, but is close for max rate
End of explanation
<END_TASK>