I built a PyMC model that predicts and successfully fits a 2D histogram to an observed 2D histogram with a Poisson likelihood using the NUTS sampler:
z_obs = pm.Poisson('z_obs', mu=hist2d_model_mean, observed=hist2d_data)
I am able to use the resulting arviz.InferenceData make trace plots, pair plots, etc. of the posterior distributions of my model parameters. However, I would now like to plot the log-likelihood of random draws of sets of model parameters from the posterior. Is the log-likelihood automatically stored in the arviz.InferenceData or easily computable with an arviz or pymc function in post-processing?
This seems like a really basic thing to do but I’m not immediately seeing how with pymc/arviz – thanks for any help!
Hmm… according to Model comparison — PyMC 5.3.1 documentation that function can only be called within a “with model” context. Whereas in my case, since my model takes a long time to run and since I have many such models to run, I saved it as a file via trace.to_netcdf('filename.nc')
. When I read it back in using arviz via trace = az.from_netcdf('filename.nc')
and then try pm.compute_log_likelihood(trace)
I get this error:
File /opt/anaconda3/envs/pymc/lib/python3.11/site-packages/pymc/stats/log_likelihood.py:63, in compute_log_likelihood(idata, var_names, extend_inferencedata, model, sample_dims, progressbar)
40 """Compute elemwise log_likelihood of model given InferenceData with posterior group
41
42 Parameters
(...)
58
59 """
61 posterior = idata["posterior"]
---> 63 model = modelcontext(model)
65 if var_names is None:
66 observed_vars = model.observed_RVs
File /opt/anaconda3/envs/pymc/lib/python3.11/site-packages/pymc/model.py:279, in modelcontext(model)
274 model = Model.get_context(error_if_none=False)
276 if model is None:
277 # TODO: This should be a ValueError, but that breaks
278 # ArviZ (and others?), so might need a deprecation.
--> 279 raise TypeError("No model on context stack.")
280 return model
TypeError: No model on context stack.
Can that function not be used on pre-existing saved files that are re-read in as arviz.InferenceData?
You need to recreate the model, as the InferenceData does not contain information about what the model / likelihood looked like.
For future reference, the solution is to add idata_kwargs={"log_likelihood": True}
to the pm.sample call when running the model. Then the resulting trace/InferenceData will have a new log_likelihood group that will be included and saved to the netcdf file, etc. I found it here Model comparison — PyMC 5.3.1 documentation and here pymc.to_inference_data — PyMC 5.3.1 documentation