<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 ABC calibration of $I_\text{Na}$ in standardised model to unified dataset.
Step1:  Initial set-up
Load experiments used for unified dataset calibration
Step2:  Combine model and experiments to produce
Step3:  Set up prior ranges for each parameter in the model.
See the modelfile for further information on specific parameters. Prepending `log_' has the effect of setting the parameter in log space.
Step4:  Run ABC calibration
Step5:  Results analysis
<ASSISTANT_TASK:>
Python Code:
import os, tempfile
import logging
import matplotlib as mpl
import matplotlib.pyplot as plt
import seaborn as sns
import numpy as np
from ionchannelABC import theoretical_population_size
from ionchannelABC import IonChannelDistance, EfficientMultivariateNormalTransition, IonChannelAcceptor
from ionchannelABC.experiment import setup
from ionchannelABC.visualization import plot_sim_results, plot_kde_matrix_custom
import myokit
from pyabc import Distribution, RV, History, ABCSMC
from pyabc.epsilon import MedianEpsilon
from pyabc.sampler import MulticoreEvalParallelSampler, SingleCoreSampler
from pyabc.populationstrategy import ConstantPopulationSize
Explanation: ABC calibration of $I_\text{Na}$ in standardised model to unified dataset.
End of explanation
from experiments.ina_sakakibara import (sakakibara_act_nyg_adjust,
                                        sakakibara_inact_nyg_adjust,
                                        sakakibara_inact_kin_nyg_adjust,
                                        sakakibara_rec_nyg_adjust)
from experiments.ina_schneider import schneider_taum_nyg_adjust
modelfile = 'models/standardised_ina.mmt'
Explanation: Initial set-up
Load experiments used for unified dataset calibration:
 - Steady-state activation [Sakakibara1992]
 - Activation time constant [Schneider1994]
 - Steady-state inactivation [Sakakibara1992]
 - Inactivation time constant (fast+slow) [Sakakibara1992]
 - Recovery time constant (fast+slow) [Sakakibara1992]
End of explanation
observations, model, summary_statistics = setup(modelfile,
                                                sakakibara_act_nyg_adjust,
                                                sakakibara_inact_nyg_adjust,
                                                schneider_taum_nyg_adjust,
                                                sakakibara_inact_kin_nyg_adjust,
                                                sakakibara_rec_nyg_adjust)
assert len(observations)==len(summary_statistics(model({})))
g = plot_sim_results(modelfile,
                     sakakibara_act_nyg_adjust,
                     schneider_taum_nyg_adjust,
                     sakakibara_inact_nyg_adjust,
                     sakakibara_inact_kin_nyg_adjust,
                     sakakibara_rec_nyg_adjust)
Explanation: Combine model and experiments to produce:
 - observations dataframe
 - model function to run experiments and return traces
 - summary statistics function to accept traces
End of explanation
limits = {'log_ina.A': (0., 1.),
          'log_ina.p_1': (1., 5.),
          'ina.p_2': (1e-7, 0.2),
          'log_ina.p_3': (-3., 1.),
          'ina.p_4': (1e-7, 0.4),
          'log_ina.p_5': (-1., 3.),
          'ina.p_6': (1e-7, 0.2),
          'log_ina.p_7': (-4., 0.),
          'ina.p_8': (1e-7, 0.2)}
prior = Distribution(**{key: RV("uniform", a, b - a)
                        for key, (a,b) in limits.items()})
# Test this works correctly with set-up functions
assert len(observations) == len(summary_statistics(model(prior.rvs())))
Explanation: Set up prior ranges for each parameter in the model.
See the modelfile for further information on specific parameters. Prepending `log_' has the effect of setting the parameter in log space.
End of explanation
db_path = ("sqlite:///" + os.path.join(tempfile.gettempdir(), "standardised_ina_unified.db"))
logging.basicConfig()
abc_logger = logging.getLogger('ABC')
abc_logger.setLevel(logging.DEBUG)
eps_logger = logging.getLogger('Epsilon')
eps_logger.setLevel(logging.DEBUG)
pop_size = theoretical_population_size(2, len(limits))
print("Theoretical minimum population size is {} particles".format(pop_size))
abc = ABCSMC(models=model,
             parameter_priors=prior,
             distance_function=IonChannelDistance(
                 exp_id=list(observations.exp_id),
                 variance=list(observations.variance),
                 delta=0.05),
             population_size=ConstantPopulationSize(1000),
             summary_statistics=summary_statistics,
             transitions=EfficientMultivariateNormalTransition(),
             eps=MedianEpsilon(initial_epsilon=100),
             sampler=MulticoreEvalParallelSampler(n_procs=16),
             acceptor=IonChannelAcceptor())
obs = observations.to_dict()['y']
obs = {str(k): v for k, v in obs.items()}
abc_id = abc.new(db_path, obs)
history = abc.run(minimum_epsilon=0., max_nr_populations=100, min_acceptance_rate=0.01)
history = abc.run(minimum_epsilon=0., max_nr_populations=100, min_acceptance_rate=0.01)
Explanation: Run ABC calibration
End of explanation
history = History(db_path)
history.all_runs() # most recent is relevant
df, w = history.get_distribution( m=0)
df.describe()
sns.set_context('poster')
mpl.rcParams['font.size'] = 14
mpl.rcParams['legend.fontsize'] = 14
g = plot_sim_results(modelfile,
                     sakakibara_act_nyg_adjust,
                     schneider_taum_nyg_adjust,
                     sakakibara_inact_nyg_adjust,
                     sakakibara_inact_kin_nyg_adjust,
                     sakakibara_rec_nyg_adjust,
                     df=df, w=w)
plt.tight_layout()
m,_,_ = myokit.load(modelfile)
sns.set_context('paper')
g = plot_kde_matrix_custom(df, w, limits=limits)
plt.tight_layout()
Explanation: Results analysis
End of explanation
<END_TASK>