138 E. Klinger and J. Hasenauer
dependence on the number of posterior modes. Since we expected the type of
kernel density estimator to influence this dependency, we probed three differ-
ent estimators: a local estimator, a global (Silverman) estimator and a cross-
validated estimator. We found that for all investigated KDE types the popu-
lation size was roughly independent of the number of posterior modesnmodes
(Fig. 2 f). The population sizes for the cross-validated KDE were larger (Fig. 2 f),
as expected due to the (generally) smaller bandwidth (see Sect. 4 ).
3.3 Model Selection for Markov Jump Process Models
Markov jump process models constitute a practically relevant class of models.
For example, they can be used to describe chemical reactions [ 8 ]. In this context,
a common task is model selection for stochastic reaction kinetics models [ 23 , 27 ].
We investigated whether the adaptive population size method could be applied
to such model selection and studied the two Markov jump process modelsm 1
andm 2 for conversion of (chemical) speciesXto speciesY,
m 1 :X+Y
k 1
−→ 2 Y, m 2 :X
k 2
−→Y,considered before in [ 5 , 18 , 25 ]. The chemical reaction kinetics were simulated
with the Gillespie algorithm [ 8 ] and representative simulations are depicted in
Fig. 3 a. The distancedbetween two trajectoriess 1 =(X 1 ,Y 1 )ands 2 =(X 2 ,Y 2 )
was defined as the absolute sum of concentration differences of speciesXevalu-
ated atN= 20 time points:
d(s 1 ,s 2 )=∑N
n=1|X 1 (tn)−X 2 (tn)|,tn=n
NT, N=20.
The simulation was run fromt 0 = 0 untilT=0.1. The initial molecule numbers
wereX(t 0 ) = 40 andY(t 0 ) = 3 for every simulation. The reference reaction
rate, used for generation of the observed datasdata,wask 1 =2.1 (log 10 k 1 =
0 .32) [ 25 ]. The priors over the log-rates log 10 k 1 and log 10 k 2 were uniform priors
log 10 ki∼U([− 2 ,2]) for both ratesi=1,2. The prior over the models was
also uniformp 0 (m 1 )=p 0 (m 2 )=1/2. The proposal densityKmodel,t(m)for
modelmat generationtwith model probabilitiesptwas given byKmodel,t(m)=
pstaypt(m)+(1−pstay)pt(m′)ifpt(m)pt(m′)>0andKmodel,t(m)=pt(m)if
pt(m)pt(m′) = 0, in whichm′denotes the other modelm′ =mandpstay=0.7.
To assess the performance and reliability of the proposed scheme for adap-
tive selection of population sizes in model selection, we performed ABC-SMC
inference for the generated artificial datasdata withECV =0.05. While the
posterior probability ofm 1 was comparable to the one ofm 2 for the first gen-
erations, we observed that as the generations progressed and the acceptance
threshold decreased, the probability ofm 1 increased as well (Fig. 3 b). In the last
generations, the posterior probability ofm 1 was close to one, resulting in the
selection of the true model. Accordingly, the parameter posterior distribution
p(log 10 k 1 |m=m 1 ) contained the parameter used to generate the artificial data