Hello, I want to migrate a GLM model to Rolling Regression Model, but I’m encountering some issues.
I don’t know how to keep all factors non-negative in MvNormal. My current approach is to use pm.math.exp
. This actually changes the distribution of the sample, and I’m not sure if it’s a good idea. At the same time, pm.Truncated
is not implemented for multivariate distributions, and MvHalfNormal
seems to be unavailable.
Here is my current code, and I’m hoping someone could help me.
signal_num = data.shape[1]-1
COORDS = {'coeffs':data.iloc[:,:-1].columns.tolist(),'time':data.index.values}
model = pm.Model(coords = COORDS)
with model:
sigma = pm.HalfNormal("sigma", sigma=0.1)
X = pm.MutableData("X",data.iloc[:,:-1].values,dims=('time','coeffs'))
beta_sigma = pm.Exponential("beta_sigma", 10.0,dims='coeffs')
beta = pm.MvGaussianRandomWalk("beta",
mu=np.zeros(signal_num),
dims=("time","coeffs"),
cov=pytensor.tensor.diag(beta_sigma),
init_dist=pm.MvNormal.dist(0.1, 0.01*np.eye(signal_num)))
Y = pm.MutableData("Y",data['Y'].values)
mu = pm.math.sum(X*pm.math.exp(beta),axis=1)
Y_obs = pm.Normal("Y_obs", mu=mu, sigma=sigma, observed=Y)
idata = jax.sample_numpyro_nuts(chains=1,draws=1000,tune=1000,postprocessing_backend="cpu",)
My original GLM model is
COORDS = {'coeffs':data.iloc[:,:-1].columns.tolist()}
model = pm.Model(coords = COORDS)
with model:
sigma = pm.HalfNormal("sigma", sigma=0.1)
X = pm.MutableData("X",data.iloc[:,:-1].values)
Y = pm.MutableData("Y",data['Y'].values)
beta = pm.TruncatedNormal(name = 'beta',mu=0.1,sigma=0.1,lower=0,dims="coeffs")
mu = pm.math.dot(X,beta)
Y_obs = pm.Normal("Y_obs", mu=mu, sigma=sigma, observed=Y)
idata = jax.sample_numpyro_nuts(chains=4,draws=10000,tune=10000)