Understanding multilevel model coefficients

We can follow your model code to see what we should expect.

The sum you calculate is for actor=0, block=0, and treatment=0. However, it’s not clear if that combination corresponds to any of your observations (and if it does, it’s not clear which observation(s) it corresponds to). If we want to see the p calculated for the first observation (i.e., d.pulled_left.values[0]), then we need to figure out what actor/block/treatment is associated with that first observation. The combination should be actor_[0] and block_[0] and treatment_[0]. Thus, the sum should be:

a[actor_[0]] + g[block_[0]] + b[treatment_[0]]

Your model immediately transforms this sum using pm.math.invlogit(). So p[0] should be

p[0] = 1 / (1 + exp(a[actor_[0]] + g[block_[0]] + b[treatment_[0]]))

Or, you can just calculate the values for all observations all at once:

p = 1 / (1 + exp(a[actor_] + g[block_] + b[treatment_]))

and check the value associated with the first observation, p[0].