Best way I would say is to use pm.Data, which is design to handle prediction better:
training_set = df_res[~df_res.index.isin(t_pred_idx)]
test_set = df_res.loc[t_pred_idx]
with pm.Model() as model_logistic:
predictor = pm.Data('predictor', training_set.t.values)
observed = pm.Data('observed', training_set.y.values)
sigma = pm.HalfNormal('sigma', sigma=5) # Half Normal are already bounded
alpha = pm.Bound(pm.Normal, lower=0.0)('alpha', mu=1, sigma=1,shape=len(np.unique(indexer))) # pm.Normal('alpha',mu=1, sigma=1) # 1,1 # make boundary: https://docs.pymc.io/api/bounds.html
beta = pm.Bound(pm.Normal, lower=0.0)('beta', mu=1, sigma=1,shape=len(np.unique(indexer))) #pm.Normal('beta',mu=1, sigma=1) # 1,1
t0 = pm.Bound(pm.Normal, lower=0.0)('t0', mu=10, sigma=10,shape=len(np.unique(indexer))) # pm.Normal('t0',mu=10, sigma=10) #normal(10,10);
mu = alpha / (1 + pm.math.exp(-(beta*(predictor-t0))))
# do indexing here instead
Y_obs = pm.Normal('Y_obs', mu=mu[indexr], sigma=sigma, observed=observed)
The nice thing here is that, by pushing the indexing to observed, for prediction you can just use mu
:
with model_logistic:
pm.set_data({'predictor': test_set.t.values})
Y_predict = pm.Normal('Y_predict', mu=mu, sigma=sigma, shape=...)
As for numerical integration, you can take a look at Custom theano Op to do numerical integration