What am I missing in this simple model?

I am working through Richard McElreath’s Statistically Rethinking.

On week 5 exercise I think I got pretty much the right answer, it certainly looks OK versus the given solution, but for some reason my model does not converge and comes up with lots of horrible warnings, whereas the given solution behaves perfectly.

After some hours poring over this, I can’t see why this would be from the two bits of code…

Given answer:

library(rethinking)
data(Wines2012)
d <- Wines2012
dat_list2b <- list(
S = standardize(d$score),
wid = d$wine.amer + 1L,
jid = d$judge.amer + 1L,
fid = ifelse(d$flight=="red",1L,2L)
)

m2b <- ulam(
alist(
S ~ dnorm( mu , sigma ),
mu <- w[wid] + j[jid] + f[fid],
w[wid] ~ dnorm( 0 , 0.5 ),
j[wid] ~ dnorm( 0 , 0.5 ),
f[wid] ~ dnorm( 0 , 0.5 ),
sigma ~ dexp(1)
), data=dat_list2b , chains=4 , cores=4 )
precis( m2b,2 )

My answer:

dat_list2 <- list(scores_s = standardize(d$score), judge_am = as.integer(d$judge.amer)+1L,
                  wine_am = as.integer(d$wine.amer)+1L, rw = as.integer(d$flight))

flist <- alist(score_s ~ dnorm(mu, sigma), 
mu <- J[judge_am] + W[wine_am] + F[rw],
J[judge_am] ~ dnorm(0, 0.5),
W[wine_am] ~ dnorm(0,0.2),
F[rw] ~ dnorm(0,0.2),
sigma ~ dexp(1))

m2 <- ulam(flist, data = dat_list2, chains = 4, cores = 4)

precis(m2, 2)

which results in horror-show below, any ideas why??

Warning messages:
1: In system(paste(CXX, ARGS), ignore.stdout = TRUE, ignore.stderr = TRUE) :
‘C:/rtools40/usr/mingw_/bin/g++’ not found
2: There were 112 divergent transitions after warmup.
to find out why this is a problem and how to eliminate them.
3: There were 1483 transitions after warmup that exceeded the maximum treedepth. Increase max_treedepth above 10.
4: There were 4 chains where the estimated Bayesian Fraction of Missing Information was low.
5: Examine the pairs() plot to diagnose sampling problems

6: The largest R-hat is 3.26, indicating chains have not mixed.
Running the chains for more iterations may help.
7: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
Running the chains for more iterations may help.
8: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
Running the chains for more iterations may help.

Hi Paul, I’m not the person to look to for actual coding errors as I use brms not Stan directly, but did you try running the model with the priors on w and f being exactly the same as Richard’s? One possibility could be that the priors you set are too constrained and so the model can’t find a nice jointly plausible set of parameters over J, W, F, and sigma? At least worth a quick try changing them from 0.2 to 0.5.

1 Like

Hi @JimBob - thanks for your response - I did try that and sadly it doesn’t solve the problem…

It looks you’re not using the same indexing as the answer. The answer has w/f/j all indexed by wid:

w[wid] ~ dnorm( 0 , 0.5 ),
j[wid] ~ dnorm( 0 , 0.5 ),
f[wid] ~ dnorm( 0 , 0.5 ),

But in your answer, they’re indexed separately:

J[judge_am] ~ dnorm(0, 0.5),
W[wine_am] ~ dnorm(0,0.2),
F[rw] ~ dnorm(0,0.2),

So I’d guess you would have to change the indexes all to wine_am

Oh yeah. Oops. Good spot @andrjohns! Guess I had just been looking at it too long. But I ran the code and the same error came up. So now I’m totally confused.
(I reran the model answer code again, and again it ran fine.) Slightly tearing my hair out…

dat_list2 <- list(scores_s = standardize(d$score), judge_am = as.integer(d$judge.amer)+1L,
                  wine_am = as.integer(d$wine.amer)+1L, rw = as.integer(d$flight))

flist <- alist(score_s ~ dnorm(mu, sigma), 
mu <- J[judge_am] + W[wine_am] + F[rw],
J[wine_am] ~ dnorm(0, 0.5),
W[wine_am] ~ dnorm(0,0.5),
F[wine_am] ~ dnorm(0,0.5),
sigma ~ dexp(1))

m2 <- ulam(flist, data = dat_list2, chains = 4, cores = 4)

Warning messages:
1: In system(paste(CXX, ARGS), ignore.stdout = TRUE, ignore.stderr = TRUE) :
‘C:/rtools40/usr/mingw_/bin/g++’ not found
2: There were 60 divergent transitions after warmup.
to find out why this is a problem and how to eliminate them.
3: There were 1543 transitions after warmup that exceeded the maximum treedepth. Increase max_treedepth above 10. See


4: There were 4 chains where the estimated Bayesian Fraction of Missing Information was low.
5: Examine the pairs() plot to diagnose sampling problems

6: The largest R-hat is 3.49, indicating chains have not mixed.
Running the chains for more iterations may help.
7: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
Running the chains for more iterations may help.
8: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
Running the chains for more iterations may help.

Try restarting R between runs. Bit drastic, but it will make sure that ulam isn’t using a cached version of the old model

Yes, I tried that too, but problem persists. Also I thought that RStan does recompiling so as not to re-use previously run models?