def compute_predictive_ll(S_test, S_train,
                          true_model=None,
                          bfgs_model=None,
                          sgd_models=None,
                          gibbs_samples=None,
                          vb_models=None,
                          svi_models=None):
    plls = {}

    # Compute homogeneous pred ll
    T = S_train.shape[0]
    T_test = S_test.shape[0]
    lam_homog = S_train.sum(axis=0) / float(T)
    plls['homog']  = 0
    plls['homog'] += -gammaln(S_test+1).sum()
    plls['homog'] += (-lam_homog * T_test).sum()
    plls['homog'] += (S_test.sum(axis=0) * np.log(lam_homog)).sum()

    if true_model is not None:
        plls['true'] = true_model.heldout_log_likelihood(S_test)

    if bfgs_model is not None:
        assert isinstance(bfgs_model, DiscreteTimeStandardHawkesModel)
        plls['bfgs'] = bfgs_model.heldout_log_likelihood(S_test)

    if sgd_models is not None:
        assert isinstance(sgd_models, list)

        plls['sgd'] = np.zeros(len(sgd_models))
        for i,sgd_model in enumerate(sgd_models):
            plls['sgd'] = sgd_model.heldout_log_likelihood(S_test)

    if gibbs_samples is not None:
        print "Computing pred ll for Gibbs"
        # Compute log(E[pred likelihood]) on second half of samplese
        offset       = len(gibbs_samples) // 2
        # Preconvolve with the Gibbs model's basis
        # F_test = gibbs_samples[0].basis.convolve_with_basis(S_test)
        S_ct, C_ct, T = convert_discrete_to_continuous(S_test, 0.02)

        plls['gibbs'] = []
        for s in gibbs_samples[offset:]:
            # plls['gibbs'].append(s.heldout_log_likelihood(S_test, F=F_test))
            plls['gibbs'].append(s.heldout_log_likelihood(S_ct, C_ct, T))

        # Convert to numpy array
        plls['gibbs'] = np.array(plls['gibbs'])

    import ipdb; ipdb.set_trace()

    if vb_models is not None:
        print "Computing pred ll for VB"
        # Compute predictive likelihood over samples from VB model
        N_models  = len(vb_models)
        N_samples = 100
        # Preconvolve with the VB model's basis
        F_test = vb_models[0].basis.convolve_with_basis(S_test)

        vb_plls = np.zeros((N_models, N_samples))
        for i, vb_model in enumerate(vb_models):
            for j in xrange(N_samples):
                vb_model.resample_from_mf()
                vb_plls[i,j] = vb_model.heldout_log_likelihood(S_test, F=F_test)

        # Compute the log of the average predicted likelihood
        plls['vb'] = -np.log(N_samples) + logsumexp(vb_plls, axis=1)

    if svi_models is not None:
        print "Computing predictive likelihood for SVI models"
        # Compute predictive likelihood over samples from VB model
        N_models  = len(svi_models)
        N_samples = 1
        # Preconvolve with the VB model's basis
        F_test = svi_models[0].basis.convolve_with_basis(S_test)

        svi_plls = np.zeros((N_models, N_samples))
        for i, svi_model in enumerate(svi_models):
            # print "Computing pred ll for SVI iteration ", i
            for j in xrange(N_samples):
                svi_model.resample_from_mf()
                svi_plls[i,j] = svi_model.heldout_log_likelihood(S_test, F=F_test)

        plls['svi'] = -np.log(N_samples) + logsumexp(svi_plls, axis=1)

    return plls
Compute the predictive log likelihood
dataset/ETHPy150Open slinderman/pyhawkes/experiments/chalearn_small_comparison.py/compute_predictive_ll