SMC Sampling

[hanzy1110]

Hi @IvanYashchuk! I recently encountered your library and wanted to thank you for the implementation its come quite handy for a problem we are tackling here at my Uni. However I'm having the following problems:

• I've read in this post about the Lokta-Volterra model that the SMC Sampler works best for the posterior distribution an ODE inverse problem has. Problem is, during my implementation when I try to use the pm.sample_smc method I get the following error (when calling pm.sample_smc(55):
  TypeError: ('Bad input argument to theano function with name "/home/user/anaconda3/envs/fenicsproject/lib/python3.9/site-packages/pymc3/distributions/distribution.py:854" at index 0 (0-based). ', 'Non-unit value on shape on a broadcastable dimension.', (55,), (True,))
  As I understand this a problem with how the output of the create_fenics_theano_op gets interpreted by the sampler. Is there a way to correct this?
• Another problem I'm facing is that when using multiple cores during sampling after sampling is complete the program just exits without any warning and no further commands are carried out. Is this a known issue?
  Thanks for your attention!!

[IvanYashchuk]

Hi @hanzy1110!
I haven't tried the SMC sampler and don't know how it works. Thanks for reporting that it doesn't work! Possibly some changes are needed to create_fenics_theano_op so that broadcasting is supported. I don't have an immediate solution to fix that.

By multiple cores do you mean running the PyMC3 sampler with cores=N or running with MPI via mpirun -n N python fenics_script.py? I guess you're talking about cores=N, I've noticed that sometimes the program hangs if N > 1 used. I didn't investigate further and haven't fixed that. What this parameter does is simply duplicating the sampler with a different random seed, so a workaround for this is to run the program multiple times (sequentially or using GNU parallel) with different random seeds save the different chains and then load with ArviZ for example to calculate the diagnostics that require multiple chains.

[hanzy1110]

Thanks for your answer! I will try your solution for the multiple cores problem but I lack enough knowledge about the innards of pymc and fenics to edit create_fenics_theano_op. Anything I should look for to start tinkering with it?