def expected_hmm_logprob(pi_0, trans_matrix, stats):
    """
    :param pi_0:          initial distribution
    :param trans_matrix:  transition matrix
    :param stats:         tuple (E[z_t], \sum_t E[z_t z_{t+1}.T])

    :return:  E_{q(z)} [ log p(z) ]
    """
    E_z, sum_E_ztztp1T, _ = stats
    T, K = E_z.shape
    assert sum_E_ztztp1T.shape == (K, K)

    out = 0
    out += np.dot(E_z[0], np.log(pi_0))
    out += np.sum(sum_E_ztztp1T * np.log(trans_matrix))
    return out