Problem with Deterministic variable

Hello, I am very new in bayesian statistics and pymc3.
I try with this very simple model

with pm.Model() as modelo_complacencia:
μ = pm.Uniform(‘μ’ , lower = 1, upper= 15)
σ = pm.Uniform(‘σ’, lower = 1, upper= 15)
P0 = pm.Normal(‘P0’, mu=μ, sd=σ, observed=data)
IPVI = pm.Deterministic(‘PVI’, get_pvi( P0))
trace_g = pm.sample(10)
az.plot_trace(trace_g)

with

def get_pvi(x):
return 1.5497/np.log10(45.22/x)

data=np.array([ 9.4634, 9.8030, 9.0072, 9.4574, 9.8042,9.0752, 9.7246, 9.9534, 9.2020, 9.8305, 9.9323, 9.2103, 9.8232, 9.7856, 9.0838, 9.9306, 9.7179, 8.9507, 9.7838, 9.4283, 8.8593, 9.8523, 9.5175, 8.9576, 9.9520, 9.4569, 9.0570, 10.1093, 9.5680, 9.1668, 10.0123, 9.4289, 9.1801, 10.1137, 9.5103, 9.2803, 10.2759, 9.4823, 9.3829, 10.3293, 9.6556, 9.6288, 10.4406, 9.6630, 9.6831, 10.5126, 9.6760, 9.7955, 10.4225, 9.5551, 9.6265, 10.2659, 9.3821, 9.7273, 10.1402, 9.3808, 9.6738, 10.4288, 9.4218, 10.0775, 10.3577, 9.6004, 10.1041, 10.2955, 9.5451, 10.2484, 10.3395, 9.5534, 10.2460, 10.1394, 9.4752, 10.3348, 10.2022, 9.5235, 10.4362, 10.1432, 9.4596, 10.5102,
10.1962, 9.7203, 10.7773, 10.3343, 9.9001, 10.8562, 10.3759, 10.0302, 11.0552, 10.2961, 9.8981, 10.8033, 10.0477, 9.7508, 10.6555, 9.9330, 9.7585, 10.5959 ])

and this message appears

anaconda3/lib/python3.8/site-packages/arviz/stats/density_utils.py:770: UserWarning: Something failed when estimating the bandwidth. Please check your data warnings.warn(“Something failed when estimating the bandwidth. Please check your data”)

Could yo please help me to understand and solve this problem?

Thanks in advance

Francisco

IPVI is always the same because it’s a deterministic function of observed data. Arviz is trying to calculate bin edges for a histogram somewhere but since all samples of IPVI are the same, the left bin edge and right bin edge are identical, so the bandwidth (measured as left minus right) is zero. In short, declaring IPVI as a deterministic doesn’t really make sense - you can compute its value outside of PyMC3. Otherwise, if P0 is not observed, then this would be a reasonable thing to do.

Chris Krapu,

Thank you very much for yor reply!!! IPVI is a function of P0, as P0 is a random variable, IPVI should be a random variable also… How can I define it?

I could define the model in this way

with pm.Model() as modelo_complacencia:
μ = pm.Uniform(‘μ’ , lower = 1, upper= 15)
σ = pm.Uniform(‘σ’, lower = 1, upper= 15)
P0 = pm.Normal(‘P0’, mu=μ, sd=σ, observed=data)
IPVI = pm.Deterministic(‘IPVI’, get_pvi(μ))
trace_g = pm.sample(10)
az.plot_trace(trace_g)

But I think in this model IPVI is a function of the mean value of the prior parameter, Is that OK? I would like to IPVI as a function of a posteriori mean parameter

Yes! You’ve got it now. That should work better and it also makes sense from a modeling perspective.

It will be a function of the p

It will be a function of the posterior mu value, since during pm.sample you are getting samples of the posterior mu

If your get_pvi is a vectorized function, you can do this much faster by calling get_pvi(trace['mu']) after your sampling is over.

Ricardo ! Thanl you very much

Dear Chris, Thank you!!!