def resample_C(self, augmented_data_list):
        """
        Resample the observation vectors. Since the emission noise is diagonal,
        we can resample the rows of C independently
        :return:
        """
        # Get the prior
        prior_precision = 1./self.sigma_C * np.eye(self.n)
        prior_mean = np.zeros(self.n)
        prior_mean_dot_precision = prior_mean.dot(prior_precision)

        # Get the sufficient statistics from the likelihood
        lkhd_precision = np.zeros((self.K-1, self.n, self.n))
        lkhd_mean_dot_precision = np.zeros((self.K-1, self.n))

        for data in augmented_data_list:
            # Compute the residual activation from other components
            # TODO: Fix this hack
            if "z" in data:
                z = data["z"]
            elif "states" in data:
                z = data["states"].stateseq
            else:
                raise Exception("Could not find latent states in augmented data!")

            # Get the observed mean and variance
            omega = data["omega"]
            kappa = data["kappa"]
            prec_obs = omega
            mu_obs = kappa / omega - self.mu[None, :]
            mu_dot_prec_obs = omega * mu_obs

            # Update the sufficient statistics for each neuron
            for k in xrange(self.K-1):
                lkhd_precision[k,:,:] += (z * prec_obs[:,k][:,None]).T.dot(z)
                lkhd_mean_dot_precision[k,:] += \
                    (mu_dot_prec_obs[:,k]).T.dot(z)

        # Sample each row of C
        for k in xrange(self.K-1):
            post_prec = prior_precision + lkhd_precision[k,:,:]
            post_cov  = np.linalg.inv(post_prec)
            post_mu   =  (prior_mean_dot_precision +
                          lkhd_mean_dot_precision[k,:]).dot(post_cov)
            post_mu   = post_mu.ravel()

            self.C[k,:] = np.random.multivariate_normal(post_mu, post_cov)
dataset/ETHPy150Open HIPS/pgmult/pgmult/distributions.py/PGMultinomialRegression.resample_C/original