Hi I was working with some hierarchical bayesian model that looks like the partially pooled hierarchical model in the eight schools example. However, I am training the model to predict different trajectories of a feature that come from several different subjects.
In my code, my feature (e.g., time) is in X, its shape is (100, 50) where 100 is the number of time samples and 50 corresponds to 50 different subjects. Also, y has the target output. This code would work when the number of time samples are equivalent between subjects. However, if I have different number of time samples per subject how would I set up the pyro.plate?
I can probably zero pad X but unless I specify the exact length of the time series, it might learn incorrect model.
Can someone suggest a solution to this problem?
kk = np.linspace(0, 10, num=100)
samp_t = np.repeat(kk,50).reshape((100,50),order='F')
samp_traj = np.zeros((100,50))
for i in range(50):
traj_types = bool(np.random.rand() > .5)
if traj_type == 1:
samp_traj[:,i] = 0.2 * samp_t[:,i]
else:
samp_traj[:,i] = 5 * samp_t[:,i]
X = torch.tensor(samp_t[:, :50], dtype=torch.float)
y = torch.tensor(samp_traj[:, :50], dtype=torch.float)
# %% Training
@config_enumerate
def model(X, y):
mu = pyro.sample("mu", dist.Normal(0., 5.))
sigma = pyro.sample("sigma", dist.Uniform(0., 10.))
number_of_traj = y.size()[1]
with pyro.plate('trajs', number_of_traj):
k = pyro.sample("k", dist.Normal(mu, sigma))
mean = k * X
bigsigma = pyro.sample("bigsigma", dist.Uniform(0., 10.))
data_length = y.size()[0]
with pyro.plate("data", data_length):
pyro.sample("obs", dist.Normal(mean, bigsigma), obs=y)
# %% Using SVI
pyro.set_rng_seed(1524)
guide = AutoDiagonalNormal(poutine.block(model, expose=[
'mu', 'sigma', 'k', 'bigsigma']))
svi = SVI(model,
guide,
optim.Adam({"lr": .01}),
loss=TraceEnum_ELBO())
pyro.clear_param_store()
num_iters = 2000
losses = []
j = -1
for i in range(num_iters):
elbo = svi.step(X, y)
losses.append(elbo)
print(elbo)
EDIT: I have tried to change my model like so:
@config_enumerate
def model(X, y, series_lengths):
mu = pyro.sample("mu", dist.Normal(0., 5.))
sigma = pyro.sample("sigma", dist.Uniform(0., 10.))
number_of_traj = y.size()[1]
k = torch.zeros(50)
bigsigma = torch.zeros(50)
mean = torch.zeros(50, 50)
for i in pyro.plate('trajs', number_of_traj):
k[i] = pyro.sample("k_{}".format(i), dist.Normal(mu, sigma))
mean[:, i] = k[i] * X[:, i]
bigsigma[i] = pyro.sample(
"bigsigma_{}".format(i), dist.Uniform(0., 10.))
data_length = series_lengths[i]
with pyro.plate("data_{}".format(i), data_length):
pyro.sample("obs_{}".format(i), dist.Normal(
mean[:data_length, i], bigsigma[i]), obs=y[:data_length, i])
# %% Using SVI
pyro.set_rng_seed(1524)
guide = AutoDiagonalNormal(poutine.block(model, expose=[
'mu', 'sigma', 'k', 'bigsigma'])) # AutoDelta(model)
svi = SVI(model,
guide,
optim.Adam({"lr": .01}),
loss=TraceEnum_ELBO())
series_lengths = np.random.randint(10, 50, 50)
pyro.clear_param_store()
num_iters = 200
losses = []
j = -1
for i in range(num_iters):
elbo = svi.step(X, y, series_lengths)
losses.append(elbo)
print(elbo)
where I send in lengths of different trajectories using an argument called series_length. But it never converges. I am sure what I have done is wrong, could someone suggest a solution please. Should I be using pyro.markov instead of outermost pyro.plate?
(but maybe there are some references showing how to tackle that issue).
