Two kinds of data in one model

Questions
2

A few thoughts:

  • You can’t use Numpy math operations on PyMC3 variables in place of Theano / PyMC3 operations (the latter are equivalent - PyMC3 just wraps the Theano code)

  • For readability, you could replace the distribution for counts_obs as

gaussian_sum = pm.math.sum(amp*pm.math.exp(-(bins-mu)**2 / sig**2)))
counts_obs = pm.Normal('counts', mu=gaussian_sum, sigma=np.sqrt(counts), observed=counts)

  • Using uniform priors is usually not a good idea - I suggest replacing all of those pm.Uniform variables with normal priors, using the distance between the uniform bounds as the standard deviation for the normal.

  • Using such a small standard deviation means that a tiny change in the parameter values leads to a huge jump in the log-likelihood which is probably why the sampler isn’t happy. Can you just work in a different set of units, i.e. if you are using joules or eV, just work in microjoules or meV? That way, you can use a sigma of 1.