def sample_W(self, nSamples, raw_samples=False):
assert isinstance(self.kern, kern.Linear)<EOL>from ..util.linalg import pdinv<EOL>N, D = self.Y.shape<EOL>Q = self.X.shape[<NUM_LIT:1>]<EOL>noise_var = self.likelihood.variance.values<EOL>Xs = np.random.randn(*((nSamples,)+self.X.shape))*np.sqrt(self.X.variance.values)+self.X.mean.values<EOL>b = np.random.rand(*((nSamples,)+self.X.shape))<EOL>Xs[b>self.X.gamma.values] = <NUM_LIT:0><EOL>invcov = (Xs[:,:,:,None]*Xs[:,:,None,:]).sum(<NUM_LIT:1>)/noise_var+np.eye(Q)<EOL>cov = np.array([pdinv(invcov[s_idx])[<NUM_LIT:0>] for s_idx in range(invcov.shape[<NUM_LIT:0>])])<EOL>Ws = np.empty((nSamples, Q, D))<EOL>tmp = (np.transpose(Xs, (<NUM_LIT:0>,<NUM_LIT:2>,<NUM_LIT:1>)).reshape(nSamples*Q,N).dot(self.Y)).reshape(nSamples,Q,D)<EOL>mean = (cov[:,:,:,None]*tmp[:,None,:,:]).sum(<NUM_LIT:2>)/noise_var<EOL>zeros = np.zeros((Q,))<EOL>for s_idx in range(Xs.shape[<NUM_LIT:0>]):<EOL><INDENT>Ws[s_idx] = (np.random.multivariate_normal(mean=zeros,cov=cov[s_idx],size=(D,))).T+mean[s_idx]<EOL><DEDENT>if raw_samples:<EOL><INDENT>return Ws<EOL><DEDENT>else:<EOL><INDENT>return Ws.mean(<NUM_LIT:0>), Ws.std(<NUM_LIT:0>)<EOL><DEDENT>
Sample the loading matrix if the kernel is linear.