You can track the cov of the full rank approximation via:

```
with model_vi:
advi = pm.FullRankADVI()
tracker = pm.callbacks.Tracker(
mean=advi.approx.params[1].eval, # callable that returns mean
L_tril=advi.approx.params[0].eval # callable that returns std
)
approx = advi.fit(50000, callbacks=[tracker])
fig = plt.figure(figsize=(16, 9))
mu_ax = fig.add_subplot(221)
std_ax = fig.add_subplot(222)
hist_ax = fig.add_subplot(212)
mu_ax.plot(tracker['mean'][8000:])
mu_ax.set_title('Mean track')
std_ax.plot(tracker['L_tril'][8000:])
std_ax.set_title('Std track')
hist_ax.plot(advi.hist[8000:])
hist_ax.set_title('Negative ELBO track');
```

to recover the cov from the lower triangular matrix:

```
n = approx.ddim
def L2full(L, n):
L_tril = np.zeros((n, n))
L_tril[np.tril_indices(n)] = L
return L_tril.dot(L_tril.T)
L_inf = np.asarray(tracker['L_tril'])[~np.isfinite(advi.hist),:]
L_inf = [L2full(L, n) for L in L_inf]
```

I had a quick look, the estimation (both `mu`

and `L_tril`

) looks quite OK. I am not completely sure where the inf comes from…