# Unexpected Lack of Convergence of Model with Multicollinearity

Hi!

I’m working on Richard Mcelreath’s book Statistical Rethinking for self-study and I’m answering end of chapter exercises using PyMC3.

One question involves creating a model with highly correlated data. The idea is to predict height with the left leg and the right leg. Data was generated as follows:

``````height = np.random.normal(loc=10, scale=2, size=100)
leg_prop = np.random.uniform(low=0.4, high=0.5, size=100)
leg_left = leg_prop * height + np.random.normal(loc=0, scale=0.02, size=100)
leg_right = leg_prop * height + np.random.normal(loc=0, scale=0.02, size=100)

highly_corr_df = pd.DataFrame({
'height': height,
'leg_left': leg_left,
'leg_right': leg_right
})
``````

The first model is simple and converges:

``````with pm.Model() as highly_corr_model:
bL = pm.Normal('bL', 2, 10)
bR = pm.Normal('bR', 2, 10)
alpha = pm.Normal('alpha', 10, 100)
sigma = pm.HalfCauchy('sigma', 1)
mu = bL * highly_corr_df['leg_left'] + bR * highly_corr_df['leg_right'] + alpha
height = pm.Normal('height', mu, sigma, observed=highly_corr_df['height'])

trace_highly_corr_model = pm.sample(4000, tune=1000, start={
'bL': 0,
'bR': 0,
'alpha': 10,
'sigma': 1
})
``````

However, the question asks what happens if `bR` was only positive. Here’s the R equivalent:

``````m5.8s2 <- map2stan(
alist(
height ~ dnorm(mu, sigma),
mu <- a + bl*leg_left + br*leg_right,
a ~ dnorm(10,100),
bl ~ dnorm(2, 10),
br ~ dnorm(2, 10) & T[0,],
sigma ~ dcauchy(0, 1),
),
data=d, chains=4,
start=list(a=10, bl=0, br=0, sigma=1)
)
``````

Here’s my code that uses `Bound` to ensure that bR can only be positive:

``````with pm.Model() as highly_corr_model_truncated:
bL = pm.Normal('bL', 2, 10, testval=0)
BoundedNormal = pm.Bound(pm.Normal, lower=0.0)
bR = BoundedNormal('bR', mu=2, sd=10, testval=0.1)
alpha = pm.Normal('alpha', 10, 100, testval=10)
sigma = pm.HalfCauchy('sigma', 1, testval=1)
mu = bL * highly_corr_df['leg_left'] + bR * highly_corr_df['leg_right'] + alpha
height = pm.Normal('height', mu, sigma, observed=highly_corr_df['height'])

trace_highly_corr_model_truncated = pm.sample(4000, tune=1000)
``````

This model has convergence issues: !

Any ideas why this is happening? I don’t know if this is a bug or if I wrote code incorrectly. Any help would be appreciated! Thanks.

The implementation looks correct to me, so the unidentifiablitly might be an issue of the model itself.
Did you try to fit the same model/data in R using `rethinking`?