code: def bsmahal(a, b, n_boot=200):
      """
      Bootstraps Mahalanobis distances for Shepherd's pi correlation.
  
      Parameters
      ----------
      a : ndarray (shape=(n, 2))
          Data
      b : ndarray (shape=(n, 2))
          Data
      n_boot : int
          Number of bootstrap samples to calculate.
  
      Returns
      -------
      m : ndarray (shape=(n,))
          Mahalanobis distance for each row in a, averaged across all the
          bootstrap resamples.
      """
      n, m = b.shape
      MD = np.zeros((n, n_boot))
      nr = np.arange(n)
      xB = np.random.choice(nr, size=(n_boot, n), replace=True)
      # Bootstrap the MD
      for i in np.arange(n_boot):
          s1 = b[xB[i, :], 0]
          s2 = b[xB[i, :], 1]
          X = np.column_stack((s1, s2))
          mu = X.mean(0)
          _, R = np.linalg.qr(X - mu)
          sol = np.linalg.solve(R.T, (a - mu).T)
          MD[:, i] = np.sum(sol**2, 0) * (n - 1)
      # Average across all bootstraps
      return MD.mean(1)
  
imports: import warnings
  import numpy
  import pandas
  import pandas_flavor
  from scipy.spatial.distance import pdist, squareform
  from scipy.stats import pearsonr, spearmanr, kendalltau
  from pingouin.config import options
  from pingouin.power import power_corr
  from pingouin.multicomp import multicomp
  from pingouin.effsize import compute_esci
  from pingouin.utils import remove_na, _perm_pval, _postprocess_dataframe
  from pingouin.bayesian import bayesfactor_pearson
label_desc: Use image node_id 3 for calling a global function with example usage: bsmahal(a, b, n_boot) and returns: MD