code: def build_mmd_func(sample1, sample2):
      """Build a PyTensor function that calculates the minimum mean discrepancy (MMD) statistic."""
  
      assert sample1.shape[1] == sample2.shape[1]
  
      s1 = pt.matrix(name="s1", shape=sample1.shape)
      s2 = pt.matrix(name="s2", shape=sample2.shape)
  
      X = np.concatenate((sample1, sample2), axis=0)
      test_ell = np.median(distance.pdist(X)) / 2
  
      K = pm.gp.cov.ExpQuad(sample1.shape[1], ls=test_ell)
      Kxx = K(s1)
      Kyy = K(s2)
      Kxy = K(s1, s2)
  
      n_x, n_y = s1.shape[0], s2.shape[0]
      mmd = (
          (pt.sum(Kxx) / (n_x * (n_x - 1)))
          + (pt.sum(Kyy) / (n_y * (n_y - 1)))
          - 2 * pt.sum(Kxy) / (n_x * n_y)
      )
  
      calc_mmd = pytensor.function(inputs=[s1, s2], outputs=mmd)
      return calc_mmd
  
imports: import arviz
  import numpy
  import pytensor
  import pytensor.tensor
  import pytest
  from scipy.spatial import distance
  import pymc
label_desc: Use image node_id 1 for calling a global function with example usage: build_mmd_func(sample1, sample2) and returns: calc_mmd