# Uber's Pyro less accurate than expected on toy example

Trying to understand the pyro example here: https://pyro.ai/examples/svi_part_i.html which starts with a Beta(10,10) prior, adds 10 Bernoulli likelihood datapoints with a 6,4 split. The analytic solution gives a posterior of Beta(16,14).

The guide starts at Beta(15,15) and after 2000 steps and 11 seconds of processing I get an estimated posterior of Beta(16.10,13.76) where I expected it to lock onto the exact answer. I also tried starting with the correct answer for the guide and processing for 20,000 steps taking 110 seconds.

$$\begin{array} \ Guide & steps & Posterior \\ Beta(15,15) & 2000 & Beta(16.10,13.76) \\ Beta(15,15) & 2000 & Beta(15.94,14.10) \\ Beta(15,15) & 2000 & Beta(16.05,13.83) \\ Beta(15,15) & 20000 & Beta(16.13,13.80) \\ Beta(16,14) & 2000 & Beta(16.04,13.97) \\ Beta(16,14) & 2000 & Beta(15.99,14.02) \\ Beta(16,14) & 20000 & Beta(15.82,14.09) \\ \end{array}$$ This is a toy problem and yet there are several issues:

• Unable to identify optimal solution
• More steps don't improve results

Are my expectations too high? Why is it unable to stick with the optimal solution?

from __future__ import print_function
import math
import os
import torch
import torch.distributions.constraints as constraints
import pyro
from pyro.infer import SVI, Trace_ELBO
import pyro.distributions as dist

# this is for running the notebook in our testing framework
smoke_test = ('CI' in os.environ)
n_steps = 2 if smoke_test else 2000

# enable validation (e.g. validate parameters of distributions)
assert pyro.__version__.startswith('0.3.1')
pyro.enable_validation(True)

# clear the param store in case we're in a REPL
pyro.clear_param_store()

# create some data with 6 observed heads and 4 observed tails
data = []
for _ in range(6):
data.append(torch.tensor(1.0))
for _ in range(4):
data.append(torch.tensor(0.0))

def model(data):
# define the hyperparameters that control the beta prior
alpha0 = torch.tensor(10.0)
beta0 = torch.tensor(10.0)
# sample f from the beta prior
f = pyro.sample("latent_fairness", dist.Beta(alpha0, beta0))
# loop over the observed data
for i in range(len(data)):
# observe datapoint i using the bernoulli likelihood
pyro.sample("obs_{}".format(i), dist.Bernoulli(f), obs=data[i])

def guide(data):
# register the two variational parameters with Pyro
# - both parameters will have initial value 15.0.
# - because we invoke constraints.positive, the optimizer
# will take gradients on the unconstrained parameters
# (which are related to the constrained parameters by a log)
alpha_q = pyro.param("alpha_q", torch.tensor(15.0),
constraint=constraints.positive)
beta_q = pyro.param("beta_q", torch.tensor(15.0),
constraint=constraints.positive)
# sample latent_fairness from the distribution Beta(alpha_q, beta_q)
pyro.sample("latent_fairness", dist.Beta(alpha_q, beta_q))

# setup the optimizer
adam_params = {"lr": 0.0005, "betas": (0.90, 0.999)}

# setup the inference algorithm
svi = SVI(model, guide, optimizer, loss=Trace_ELBO())

for step in range(n_steps):
svi.step(data)
if step % 100 == 0:
print('.', end='')

# grab the learned variational parameters
alpha_q = pyro.param("alpha_q").item()
beta_q = pyro.param("beta_q").item()

# here we use some facts about the beta distribution
# compute the inferred mean of the coin's fairness
inferred_mean = alpha_q / (alpha_q + beta_q)
# compute inferred standard deviation
factor = beta_q / (alpha_q * (1.0 + alpha_q + beta_q))
inferred_std = inferred_mean * math.sqrt(factor)

print("\nbased on the data and our prior belief, the fairness " +
"of the coin is %.3f +- %.3f" % (inferred_mean, inferred_std))


I lack the reputation to create a pyro tag.

• Have you tried more iterations? – Tim Apr 19 at 15:13
• 20000 steps moved the result further from the truth than 2000 as shown in the table. – Andrew Apr 20 at 4:07

As always when using SGD, care has to be taken with the learning rate. By decaying learning rate from 1E-3 to 8E-5 over 5000 iterations ($$\gamma = 0.995$$) using the SGD optimizer with momentum 0.9, I was able to get these results from 3 runs of the algorithm: for $$\alpha$$ and $$\beta$$ respectively:
15.997169494628906 13.988293647766113

Iteration count divided by 10 on the x-axis, $$\beta_q$$ on the $$y$$ axis. As the learning rate is still decaying, more iterations would continue to improve the convergence. 