code: def multiple_mahalanobis(effect, covariance):
      """Returns the squared Mahalanobis distance for a given set of samples
  
      Parameters
      ----------
      effect: array of shape (n_features, n_samples),
              Each column represents a vector to be evaluated
      covariance: array of shape (n_features, n_features, n_samples),
                  Corresponding covariance models stacked along the last axis
  
      Returns
      -------
      sqd: array of shape (n_samples,)
           the squared distances (one per sample)
      """
      # check size
      if effect.ndim == 1:
          effect = effect[:, np.newaxis]
      if covariance.ndim == 2:
          covariance = covariance[:, :, np.newaxis]
      if effect.shape[0] != covariance.shape[0]:
          raise ValueError(
              "Inconsistant shape for effect and covariance"
          )
      if covariance.shape[0] != covariance.shape[1]:
          raise ValueError("Inconsistant shape for covariance")
  
      # transpose and make contuguous for the sake of speed
      Xt, Kt = np.ascontiguousarray(effect.T), np.ascontiguousarray(
          covariance.T
      )
  
      # compute the inverse of the covariances
      Kt = multiple_fast_inv(Kt)
  
      # derive the squared Mahalanobis distances
      sqd = np.sum(
          np.sum(Xt[:, :, np.newaxis] * Xt[:, np.newaxis] * Kt, 1), 1
      )
      return sqd
  
file_path: nipy/nipy/algorithms/statistics/utils.py
imports: from itertools import combinations
  import numpy
  from scipy.stats import norm
label_desc: Use image node_id 3 for calling a global function with example usage: multiple_mahalanobis(effect, covariance) and returns: sqd