code: def _apply_matrix_kernel(
      kernel_fct, X, Y=None, *, Cref=None, reg=1e-10
  ):
      """Apply a matrix kernel function."""
      _check_dimensions(X, Y, Cref)
      n_matrices_X, n, n = X.shape
  
      X_ = kernel_fct(X, Cref)
  
      if isinstance(Y, type(None)) or np.array_equal(X, Y):
          Y_ = X_
      else:
          Y_ = kernel_fct(Y, Cref)
  
      # calculate scalar products: K[i,j] = np.trace(X_[i]^T @ Y_[j])
      X_T = X_.transpose((0, 2, 1))
      K = np.einsum("acb,dbc->ad", X_T, Y_, optimize=True)
  
      # regularization due to numerical errors
      if np.array_equal(X_, Y_):
          K.flat[:: n_matrices_X + 1] += reg
  
      return K
  
file_path: pyRiemann/pyriemann/utils/kernel.py
imports: import numpy
  from .base import invsqrtm, logm
  from .mean import mean_riemann
label_desc: Use image node_id 5 for calling a global function with example usage: _apply_matrix_kernel(kernel_fct, X, Y) and returns: K