'projection matrix (nobs_i, nobs_i) (M in regression ?)  (JP check, guessing)
Display (3.10) from Laird, Lange, Stram (see help(Unit))
W - W X Sinv X\' W\''

def compute_P(self, Sinv):

    t = np.dot(self.W, self.X)
    self.P = (self.W - np.dot(np.dot(t, Sinv), t.T))
