def rvs(self, z, N=1, full_output=False):
        T,D = z.shape
        psis = z.dot(self.C.T) + self.mu[None, :]
        pis = np.zeros((T, self.K))
        xs = np.zeros((T, self.K))
        for t in xrange(T):
            pis[t,:] = psi_to_pi(psis[t,:])
            xs[t,:] = np.random.multinomial(N, pis[t,:])

        if full_output:
            return pis, xs
        else:
            return xs
Sample from a PG augmented multinomial distribution
dataset/ETHPy150Open HIPS/pgmult/pgmult/distributions.py/PGMultinomialRegression.rvs