code: def resid_dev(self, endog, mu, scale=1.):
          
  
          return (endog - mu) / np.sqrt(self.variance(mu)) / scale
docs: Gaussian deviance residuals
  
          Parameters
          -----------
          endog : array-like
              Endogenous response variable
          mu : array-like
              Fitted mean response variable
          scale : float, optional
              An optional argument to divide the residuals by scale. The default
              is 1.
  
          Returns
          -------
          resid_dev : array
              Deviance residuals as defined below