Hello, I have trained a Bayesian forward model using polynomial regression (with input X of shape (100,10) and output Y of shape (100,15)) that results in the unknown coefficients (theta) as posterior samples. Now, I want to formulate an inverse problem where the median value is chosen from the posterior samples of theta and the inference is to be done on the random input vector X. The problem is after getting the trace, and setting the ‘y’ to an unseen target response, the posterior samples do not change much. That is to say, they are very close to what prior distributions I choose on X; and the posterior does not change. Please refer the code and suggest what could be done. Where am I going wrong? Thanks alot.
def forward_model(X):
X = scaleX(X)
temp = alpha0_med+ pt.dot(alpha1_med.T, X.T) + pt.dot(alpha2_med.T, pt.power(X,2).T) + pt.dot(alpha3_med.T,
pt.power(X,3).T)
return temp
def scaleX(x):
min_vals = pt.as_tensor_variable(scaler_x.data_min_)
scale = pt.as_tensor_variable(scaler_x.scale_)
return (x - min_vals) * scale
with pm.Model() as model_inv:
model_inv.add_coord("n", np.arange(y_training.shape[0]), mutable=True)
model_inv.add_coord("p", np.arange(10), mutable=True)
model_inv.add_coord("f", np.arange(y_training.shape[1]), mutable=True)
y = pm.MutableData("y", y_training, dims=["n", "f"])
X_rand = []
randVars_name= []
for i in range(X_data_df.shape[1]):
if i<9:
x_i = pm.Normal('x' + '_' + str(i),mu=X_data_df.iloc[:,i].mean(),sigma=X_data_df.iloc[:,i].std())
elif i==9:
x_i = pm.InverseGamma('x' + '_' + str(i), mu=X_data_df.iloc[:, i].mean(),
sigma=2)
else:
x_i = pm.InverseGamma('x' + '_' + str(i), mu=200, sigma=10)
X_rand.append(x_i)
randVars_name.append(x_i.name)
X_rand_stack = pt.stack(X_rand).T
sigma_y_ip = pm.InverseGamma("sigma_y", mu=1, sigma=1,dims="f")
likelihood = pm.Normal('likelihood', mu=forward_model(X_rand_stack), sigma=sigma_y_ip,
observed=y)
trace3 = pm.sample(draws=2000,
tune=500,
init="adapt_diag",
chains=4,
cores=4)
'''
Model validation--> To get the predictions of X given the unseen target y
'''
with model_inv:
pm.set_data(
{
'y': y_test[1:2, :],
},
coords={"n": np.arange(y_training.shape[0],y_training.shape[0] + 1)},
)
ppc_samples = pm.sample_posterior_predictive(trace3,var_names=randVars_name,predictions=True)
ppc_samples_allVar = pm.sample_posterior_predictive(trace3,predictions=True)
features = randVars_name
ncol = 3
fig, axes = plt.subplots(nrows=6,ncols=ncol, dpi=150,figsize=(12,6))
for i in range(X_data_df.shape[1]):
if i<9:
sns.distplot(ppc_samples.predictions['x_' + str(i)].values.flatten(), kde=False, hist=True,
ax=axes[i // ncol, i % ncol],
color='y')
axes[i // ncol, i % ncol].set_xlabel(f'{features[i]}')
else:
sns.distplot(ppc_samples.predictions['x_' + str(i)].values.flatten(), kde=False, hist=True,
ax=axes[i // ncol, i % ncol],
color='g')
axes[i//ncol,i%ncol].set_xlabel(f'{features[i]}')
axes[i//ncol,i%ncol].axvline(x=ppc_samples.predictions['x_'+str(i)].mean(),c='b', label='Posterior mean ('
'computed)')
axes[i//ncol,i%ncol].axvline(x=X_test.iloc[1,i], c='r', label='Known test data')
# axes[i//ncol,i%ncol].legend(fontsize=8)
axes[i // ncol, i % ncol].set_ylabel('Occurrences',fontsize=8)
axes[i // ncol, i % ncol].tick_params(axis='x', labelsize=6)