def resample(self, data, mask=None, omega=None):
        """
        Multinomial regression is somewhat special. We have to compute the
        kappa functions for the entire dataset, not just for one column of
        the data at a time.
        """
        if not isinstance(data, list):
            assert isinstance(data, tuple) and len(data) == 2, \
                "datas must be an (x,y) tuple or a list of such tuples"
            data = [data]

        if mask is None:
            mask = [np.ones(y.shape, dtype=bool) for x, y in data]

        # Resample auxiliary variables if they are not given
        if omega is None:
            omega = self._resample_auxiliary_variables(data)

        # Make copies of parameters (for sample collection in calling methods)
        self.A = self.A.copy()
        self.b = self.b.copy()

        D = self.D_in
        for n in range(self.D_out):
            # Resample C_{n,:} given z, omega[:,n], and kappa[:,n]
            prior_Sigma = np.zeros((D + 1, D + 1))
            prior_Sigma[:D, :D] = self.sigmasq_A[n]
            prior_Sigma[D, D] = self.sigmasq_b[n]
            prior_J = np.linalg.inv(prior_Sigma)

            prior_h = prior_J.dot(np.concatenate((self.mu_A[n], [self.mu_b[n]])))

            lkhd_h = np.zeros(D + 1)
            lkhd_J = np.zeros((D + 1, D + 1))

            for d, m, o in zip(data, mask, omega):
                if isinstance(d, tuple):
                    x, y = d
                else:
                    x, y = d[:, :D], d[:, D:]
                augx = np.hstack((x, np.ones((x.shape[0], 1))))
                J = o * m
                h = self.kappa_func(y) * m

                lkhd_J += (augx * J[:, n][:, None]).T.dot(augx)
                lkhd_h += h[:, n].T.dot(augx)

            post_h = prior_h + lkhd_h
            post_J = prior_J + lkhd_J

            joint_sample = sample_gaussian(J=post_J, h=post_h)
            self.A[n, :] = joint_sample[:D]
            self.b[n] = joint_sample[D]
def resample ( self , data , mask = None , omega = None ) : if not isinstance ( data , list ) : assert isinstance ( data , tuple ) and len ( data ) == 2 , "datas must be an (x,y) tuple or a list of such tuples" data = [ data ] if mask is None : mask = [ np . ones ( y . shape , dtype = bool ) for x , y in data ] if omega is None : omega = self . _resample_auxiliary_variables ( data ) self . A = self . A . copy ( ) self . b = self . b . copy ( ) D = self . D_in for n in range ( self . D_out ) : prior_Sigma = np . zeros ( ( D + 1 , D + 1 ) ) prior_Sigma [ : D , : D ] = self . sigmasq_A [ n ] prior_Sigma [ D , D ] = self . sigmasq_b [ n ] prior_J = np . linalg . inv ( prior_Sigma ) prior_h = prior_J . dot ( np . concatenate ( ( self . mu_A [ n ] , [ self . mu_b [ n ] ] ) ) ) lkhd_h = np . zeros ( D + 1 ) lkhd_J = np . zeros ( ( D + 1 , D + 1 ) ) for d , m , o in zip ( data , mask , omega ) : if isinstance ( d , tuple ) : x , y = d else : x , y = d [ : , : D ] , d [ : , D : ] augx = np . hstack ( ( x , np . ones ( ( x . shape [ 0 ] , 1 ) ) ) ) J = o * m h = self . kappa_func ( y ) * m lkhd_J += ( augx * J [ : , n ] [ : , None ] ) . T . dot ( augx ) lkhd_h += h [ : , n ] . T . dot ( augx ) post_h = prior_h + lkhd_h post_J = prior_J + lkhd_J joint_sample = sample_gaussian ( J = post_J , h = post_h ) self . A [ n , : ] = joint_sample [ : D ] self . b [ n ] = joint_sample [ D ]
Multinomial regression is somewhat special. We have to compute the
        kappa functions for the entire dataset, not just for one column of
        the data at a time.