I ran into this piece of code:
data = np.repeat((0, 1), (3, 6))
with pm.Model() as normal_approximation:
p = pm.Uniform('p', 0, 1)
w = pm.Binomial('w', n=len(data), p=p, observed=data.sum())
mean_q = pm.find_MAP() # MAP: maximum a posteriori probability
std_q = ((1/pm.find_hessian(mean_q, vars=[p]))**0.5)[0]
mean_q['p'], std_q
I tried to search google generally and look at the documentation at what is the argument observed for pm.Binomial and what is the function find_hessian (for example what's the meaning of the vars keyword) - but I see no explanation.
can someone not only explain these questions to me, but refer me to the right source to answer such questions myself?
I'd generally recommend reading through the tutorials, especially the "API Quick Start" and "Getting Started" ones. They'll cover the basics of how PyMC3 builds a graphical model of RandomVariables as a compute graph, including what observed RVs are (i.e., the ones you have observed and therefore have data from which a likelihood can be evaluated). Unless one already has a strong background in Bayesian inference, learning through examples will likely be the most effective route.
Otherwise, keep the API at hand, and make use of Python's help function.
For example, you can get the signature of find_hessian with the latter:
> help(pm.find_hessian)
find_hessian(point, vars=None, model=None)
Returns Hessian of logp at the point passed.
Parameters
----------
model: Model (optional if in `with` context)
point: dict
vars: list
Variables for which Hessian is to be calculated.
This particular method is more of an internal utility (e.g., it gets used by NUTS sampler), but those with a math-stats or physics background (such as the person who wrote the original code) will recognize it as the matrix of second partials, and understand its relationship to the covariance matrix. In the particular code, it is evaluating the p variable at the MAP location (mean_q), and then this is converted to the standard deviation.