For the project i have been working on i am trying to implement a model which is a sum of lognormals and so i wanted to test how Pymc3 handles a general lognormal of the form
here is my implementation in Pymc
with pm.Model() as model: mval = pm.Normal('mval', mu=5.0, sigma=5.0) sval = pm.Normal('sval',mu=2.0,sigma=2.0) tval = pm.Normal('tval', mu=5.0, sigma=5.0) LN = tt.exp(-tt.log((x-tval)/mval)**2/(2*sval**2))/((x-tval)*sval*tt.sqrt(2*np.pi)) LN = tt.switch(tt.isnan(LN), 0.0, LN) y_ = pm.Normal("y", mu=LN, observed=lognorm) trace = pm.sample(1000,tune=1000)
the switch is there so i do not get NaN values when tval>x. I also took a note from the Stan people and used normal priors instead of uniform ones .“lognorm” is obtained from using the same equation with the parameters mval=2.5, sval=1.5 and tval=7.5. The number of points here is set using numpy.linspace
When i run the sampler i am surprised to find that is diverges between 60-90 times for each chain after tuning. The mean parameter values found are off from the ones i inputted for lognorm and with a relatively high sd. I also find that the number of effective samples is low (less than 10%)
Can someone help me understand why the sampler is having trouble with this lognormal? Again i am surprised by the result because a lognormal only has 3 parameters. Is there some way to define the lognormal differently so it keeps the sampler from diverging?