ValueError: Bad initial energy: inf. The model might be misspecified

Hi, I have the following model but I am getting the bad initial energy error. I checked similar questions but cannot fix the problem. Any pointers would be greatly appreciated. Basically I am trying to use the radon example with some additional variables. I have panel data.

with pm.Model() as hierarchical_model:
mu_constant = pm.Normal(‘mean_for_constants’, mu=0., sd= 1)
sigma_constant = pm.Uniform(‘sigma_for_constants’, lower = 0, upper = 1)

mu_tr_coefs = pm.Normal('mean_for_trend_coefs', mu=0., sd= 1)
sigma_tr_coefs = pm.Uniform('sigma_for_tr_coefs', lower = 0, upper = 1)

mu_var_coefs = pm.Normal('mean_for_var_coefs', mu=0., sd= 1)
sigma_var_coefs = pm.Uniform('sigma_for_var_coefs', lower = 0, upper = 1)

#intercept for geos
c = pm.Normal('constant', mu=mu_constant, sd=sigma_constant, shape=len(data.geo.unique()))
#coefs for trend
b = pm.Normal('trend_coefs', mu=mu_tr_coefs, sd=sigma_tr_coefs, shape=len(data.geo.unique()))
#coefs for var
g = pm.Normal('var_coefs', mu=mu_var_coefs, sd=sigma_var_coefs, shape=len(data.geo.unique()))
epsilon = pm.Normal('epsilon', mu=0., sd= 100)

exp_value = c[geo_ids] + b[geo_ids]*data.time.values + g[geo_ids]*data.Total_var.values

outcome_likelihood = pm.Normal('outcome_likelihood', mu = exp_value, sd = epsilon, observed = data.output)

Here is the result of model.test_point

(‘mean_for_constants’, array(-5.524108719192764))
(‘sigma_for_constants_interval__’, array(-1.3862943611198906))
(‘mean_for_trend_coefs’, array(-0.9189385332046727))
(‘sigma_for_tr_coefs_interval__’, array(-1.3862943611198906))
(‘mean_for_var_coefs’, array(-0.9189385332046727))
(‘sigma_for_var_coefs_interval__’, array(-1.3862943611198906))
(‘constant’, array(-1014.5019231128882))
(‘trend_coefs’, array(-1014.5019231128882))
(‘var_coefs’, array(-1014.5019231128882))
(‘epsilon’, array(-5.524108719192764))
(‘outcome_likelihood_missing’, array(0))
(‘outcome_likelihood’, array(-inf))

Any help is greatly appreciated.

This is what causing the issue - you should not use epsilon = pm.Normal('epsilon', mu=0., sd= 100) as sigma for Normal distribution should be larger than 0. Try epsilon = pm.HalfNormal('epsilon', sd= 100)