umbertopicchini / SAEM-ABC

example code for the paper by Picchini & Samson 2017, "Coupling stochastic EM and Approximate Bayesian Computation for parameter inference in state-space models", arXiv:1512.04831.
GNU Lesser General Public License v3.0
3 stars 5 forks source link

Error in theophylline SAEM ABC #1

Open finmod opened 6 years ago

finmod commented 6 years ago

In theophylline SAEM ABC I get this error. Everything else runs fine.

Iteration #2, now using ABC tolerance = 5.000000e-01 Undefined function or variable 'XSum'.

Error in abcsmc_filter (line 21) logweights = (logweights-max(logweights))-log(XSum(exp(logweights-max(logweights)))); % normalize weights; suggestion from page 6 of Cappe et al. "An overview of existing methods and recent advances in sequential Monte Carlo"

Error in saem_abcsmc (line 39) [xhat_selected_big,xhat_selected_small] = abcsmc_filter(model_param,yobs,initstate,numparticles,N_threshold,abc_tolerance,verbose);

Error in theophylline_run (line 94) THETAmatrix_saem = saem_abcsmc(model_param,parmask,parbase,yobs,initstate,saem_numit,warmup,fisherestim_iter,numparticles,N_threshold,abc_schedule,abc_vector,verbose);

umbertopicchini commented 6 years ago

Thank you for spotting the mistake. I haven't checked yet if it is possible to replace all "XSum" entries with the regular "sum" without affecting the quality of the results. I will look into this when I have time. Therefore, in the meantime you could: 1) the easiest is to replace all entries XSum --> sum, but then I do not know if results change. If they do not change, I will reupload a version with the necessary modification. 2) (i) if you are a Windows user look into the newly uploaded files. One of those is a Matlab wrapper to XSum.c and an efficient mex version I compiled for a Windows x64 architecture is also available (XSum.mexw64). So if you use Windows you are good to go. (ii) if you are not a Windows user you will need to compile XSum.c to obtain a mex file suitable for your architecture.

Further info on XSum are available at https://se.mathworks.com/matlabcentral/fileexchange/26800-xsum?s_tid=gn_loc_drop

finmod commented 6 years ago

For once, I am on the lucky side of the OS as I am using Windows. Solution 2 works spot on. Just add these instructions to the README.

finmod commented 6 years ago

Do you have the Lotka Volterra ODE application in any of the SAEM repositories here?

@ChrisRauckaukas

umbertopicchini commented 6 years ago

There is no LV application for SAEM-ABC. You are talking of my other package, SAEM-SL. Anyway: the LV application I considered in SAEM-SL does NOT have an ODE model, underlying the dynamics. The model implemented considers noisy observations (Y) of a Markov jump process (X), the latter taking values on the set of non-negative integers. The true underlying system (X) has exact dynamics given by the so-called "Gillespie algorithm", the latter producing EXACT realizations of the solution process of a stochastic kinetic model (it is possible to interpret two interacting populations as chemical reactions between two reactants).