In the Bayesian Cognitive Modeling, chapter 12.2 : Psychophysical functions under contamination. When I run the model2b_, it gives error: AttributeError: ‘MeanFieldGroup’ object has no attribute ‘gbij’.
start = trace2b[0]
start[‘zij’] = start[‘zij’].astype(int)
stds = approx.gbij.rmap(approx.std.eval()) # cause by this
cov = model2b.dict_to_array(stds) ** 2
Could you hlpe me to solve it ?
Do stds = approx.bij.rmap(approx.params[1].eval()) instead. See also Extracting variable means from ADVI fitted model
Thank you, peng ge! 
Hi, In the Bayesian Cognitive Modeling, chapter 12.2, I want to use this model which considered the contaminant data into my own model. But the problem is that it’s so hard to converge.
You written a sentence in the file which is ‘NOTICE: In the model following I set the prior of σσ to σα,σβ,σϕ∼Uniform(0,3)σα,σβ,σϕ∼Uniform(0,3) . Otherwise the model does not converge to the desired parameters.’.
I tried the ways only use NUTS or ADVI to initialize my model just like chapter 12.2. But it doesn’t work.
So, how can I debug my own model with this contaminant data structure?
The containment data model needs some update actually, as a start, uniform prior for sigma is not recommended (works okish in Gibbs sampling, but for nuts it is not good) try halfnormal or halfcachy with a small standard deviation.
The second point to improve is rewriting the part involving the latent discrete variable. The STAN solution for these kind of model is marginalized out the discrete variable, which leaves only continuous variables in the model. You can have a look at the STAN code comes with the book about this solution.
Thank you, but I do not understand R, 
I tried halfnormal or halfcachyin in my model, the result is similar with the uniform. And my model can converge without the containment data model.
My code is like this:
start = trace_2[0]
start['zij'] = start['zij'].astype(int)
stds = approx.bij.rmap(approx.std.eval())
cov = model_2.dict_to_array(stds) ** 2
with pm.Model() as model_2b:
b = 2
alpha = pm.HalfNormal('alpha', 0.2)
mu_4 = pm.Normal('mu_4', mu=0, sd=100)
sd_4 = pm.HalfNormal('sd_4', b)
mu_3 = pm.Normal('mu_3', mu=0, sd=100)
sd_3 = pm.HalfNormal('sd_3', b)
mu_2 = pm.Normal('mu_2', mu=0, sd=100)
sd_2 = pm.HalfNormal('sd_2', b)
mu_1 = pm.Normal('mu_1', mu=0, sd=100)
sd_1 = pm.HalfNormal('sd_1', b)
mu_0 = pm.Normal('mu_0', mu=0, sd=100)
sd_0 = pm.HalfNormal('sd_0', b)
beta4 = pm.Normal('beta4', mu_4, sd_4, shape=companiesABC)
beta3 = pm.Normal('beta3', mu_3, sd_3, shape=companiesABC)
beta2 = pm.Normal('beta2', mu_2, sd_2, shape=companiesABC)
beta1 = pm.Normal('beta1', mu_1, sd_1, shape=companiesABC)
beta = pm.Normal('beta', mu_0, sd_0, shape=companiesABC)
liner = pm.Deterministic('liner', tt.exp(beta[ip] + \
(beta1[Num_shared] * xs_year + beta2[Num_shared] * xs_char1 +\
beta3[Num_shared] * xs_char2 + beta4[Num_shared] * xs_year * xs_year)))
pi_ij = pm.Uniform('pi_ij', lower=0, upper=1, shape=len(Num_shared.get_value()))
# latent model for contamination
sigma_p = pm.HalfNormal('sigma_p', b)
mu_p = pm.Normal('mu_p', mu=0, tau=.001)
probitphi = pm.Normal('probitphi', mu=mu_p, sd=sigma_p, shape=companiesABC, testval=np.ones(companiesABC))
phii = pm.Deterministic('phii', Phi(probitphi))
zij = pm.Bernoulli('zij', p=phii[Num_shared], shape=len(Num_shared.get_value()))
beta_mu = pm.Deterministic('beta_mu', tt.switch(tt.eq(zij, 0), liner, pi_ij))
Observed = pm.Weibull("Observed", alpha=alpha, beta=beta_mu, observed=ys_faults) # 观测值
step = pm.NUTS(scaling=cov, is_cov=True)
trace_2b = pm.sample(3000, step=[step], start=start, turn =500)
And if I use uniform to replace phii , than the containment data model doesn’t work at all, the result is the same as I do not use containment data model. Is there any ways? 
This is the results:
Looking at the trace, the zij never moves, and the model is either not converge or has multi mode. These are indication that the containment model is not working. Did u plot the raw data to make sure using a containment model makes sense?
