Hello,
I am trying to apply a simple linear regression model to some data with uniform likelihood, here is my code:
def linear_model_uniform(X, y,x_sigma,y_sigma,intercept_prior,coefficient_prior):
'''
A function to define a linear model in pyro
------------Inputs--------------
X: 2D torch tensor with shape (n_samples,n_features)
y: 1D torch tensor with shape (n_samples)
x_sigma: float, standard deviation of the error for age, which is obtained from the age data model
y_sigma: float, standard deviation or covariance function of the error for the RSL, which is obtained from the RSL data model
intercept_prior: pyro distribution for the intercept coefficient
coefficient_prior: pyro distribution for the slope coefficient
'''
# Define our intercept prior
linear_combination = pyro.sample("b", intercept_prior)
#Define our coefficient prior
beta_coef = pyro.sample("a", coefficient_prior)
#generate random error for age
N = X.shape[0]
x_noise = pyro.sample('obs_xerr',dist.Normal(torch.zeros(N),x_sigma).to_event(1))
x_noisy = X[:, 0]+x_noise
#calculate mean prediction
mean = linear_combination + (x_noisy * beta_coef)
with pyro.plate("data", X.shape[0]):
# Condition the expected mean on the observed target y
y_up = mean + 2*y_sigma
y_down = mean - 2*y_sigma
pyro.sample("obs", dist.Uniform(y_down,y_up), obs=y)
from pyro.infer import MCMC, NUTS,HMC
#------Define the process mdoel---------
linear_model = linear_model_uniform
test_X = torch.tensor(t_line)[:,None] #convert X to a 2D array which is more common for pytorch models
#-------Define parameter model---------
#you can change the prior distribution here to see how it affects the model!
intercept_prior = dist.Uniform(-50., 50.)
coefficient_prior = dist.Uniform(-1e-1,0)
data_prior = dist.Uniform(torch.tensor(omega+nu_l),torch.tensor(omega+nu_u))
y_sigma = torch.tensor((nu_l+nu_u)/2)
nuts_kernel = NUTS(linear_model)
mcmc = MCMC(nuts_kernel, num_samples=1000, warmup_steps=100)
mcmc.run(test_X,torch.tensor(omega),torch.tensor(t_sigma),y_sigma,intercept_prior,coefficient_prior)
#-------get posterior samples---------
hmc_samples = {k: v.detach().cpu().numpy() for k, v in mcmc.get_samples().items()}
And I always get error message like: ValueError: Model specification seems incorrect - cannot find valid initial params.
I just wondering am I doing this in the right way, and what if I want to use other distribution as my likelihood like lognormal?
Thanks!!