Model with multiple likelihoods - how do they contribute to overall log-likelihood during sampling?

[jonsedar]

Probably a stupid question, and this isn’t quite my situation but as close as I think I can get while keeping it simple.

Lets say I have a model where I think a common set of X features correlate to 2 observeds, and I create 2 likelihoods, e.g.:

coords = dict(x_jv = dfx.columns.drop(['y1', 'y2']).values)
coords_mutable = dict(oid = dfx.index.values)

with pm.Model(coords=coords, coords_mutable=coords_mutable) as mdl:
    y1 = pm.MutableData('y1', dfx['y1'].values, dims='oid')
    y2 = pm.MutableData('y2', dfx['y2'].values, dims='oid')
    x = pm.MutableData('x', dfx.drop(['y1', 'y2'], axis=1).values, dims=('oid', 'x_jv'))
    
    # 1. Create linear models
    b1 = pm.Normal('b1', mu=0.0, sigma=1.0, dims='x_jv')
    mu1 = pt.dot(b1, x.T)
    sigma1 = pm.InverseGamma('sigma1', alpha=5.0, beta=4.0)

    b2 = pm.Normal('b2', mu=0.0, sigma=1.0, dims='x_jv')
    mu2 = pt.dot(b2, x.T)
    sigma2 = pm.InverseGamma('sigma2', alpha=5.0, beta=4.0)

    # 2. Condition using observed
    _ = pm.LogNormal('yhat1', mu=mu1, sigma=sigma1, observed=y1, dims='oid')
    _ = pm.LogNormal('yhat2', mu=mu2, sigma=sigma2, observed=y2, dims='oid')

How are these two independent likelihoods treated during sampling? Are their log-likelihoods summed together?

[ricardoV94]

@jonsedar

Are their log-likelihoods summed together?

Yes.

The two following models are identical:

with pm.Model() as m1:
  mus = pm.Normal("x", shape=(2,))
  y = pm.Normal("y", mu=mus, observed=[0, 1])

with pm.Model() as m2:
  mus = pm.Normal("x", shape=(2,))
  y1 = pm.Normal("y1", mu=mus[0], observed=[0])
  y2 = pm.Normal("y2", mu=mus[1], observed=[1])

[jonsedar]

Thanks @ricardoV94, told you it was stupid question

Where I was hoping to lead to is around how the pointwise log-likelihood is stored, and now I see that they’re simply stored with increasing dimensions as needed

e.g.

with pm.Model() as m1:
  mu = pm.Normal("mu", shape=(2,))
  y = pm.Normal("y", mu=mu, observed=[0, 1])
  idata_m1 = pm.sample(idata_kwargs={"log_likelihood": True})

idata_m1.log_likelihood

with pm.Model() as m2:
  mu = pm.Normal("mu", shape=(2,))
  y1 = pm.Normal("y1", mu=mu[0], observed=[0])
  y2 = pm.Normal("y2", mu=mu[1], observed=[1])

  idata_m2 = pm.sample(idata_kwargs={"log_likelihood": True})

idata_m2.log_likelihood