code: def _compute(X, Y, settings, init_vars):
      """Compute HDMR"""
      N, d, maxorder, maxiter, m, K, R, alpha, lambdax = settings
      Em, idx, SA, RT, Y_em, Y_id, m1, m2, m3, j1, j2, j3 = init_vars
  
      # DYNAMIC PART: Bootstrap for confidence intervals
      for k in range(K):
          tic = time.time()  # Start timer
  
          # Extract the "right" Y values
          Y_id[:] = Y[idx[:, k]]
  
          # Compute the variance of the model output
          SA["V_Y"][k, 0] = np.var(Y_id)
  
          # Mean of the output
          Em["f0"][k] = np.sum(Y_id) / R
  
          # Compute residuals
          Y_res = Y_id - Em["f0"][k]
  
          # 1st order component functions: ind/backfitting
          Y_em[:, j1], Y_res, Em["C1"][:, :, k] = _first_order(
              Em["B1"][idx[:, k], :, :],
              Y_res,
              Em["C1"][:, :, k],
              R,
              Em["n1"],
              m1,
              maxiter,
              lambdax,
          )
  
          # 2nd order component functions: individual
          if maxorder > 1:
              Y_em[:, j2], Y_res, Em["C2"][:, :, k] = _second_order(
                  Em["B2"][idx[:, k], :, :],
                  Y_res,
                  Em["C2"][:, :, k],
                  R,
                  Em["n2"],
                  m2,
                  lambdax,
              )
  
          # 3rd order component functions: individual
          if maxorder == 3:
              Y_em[:, j3], Em["C3"][:, :, k] = _third_order(
                  Em["B3"][idx[:, k], :, :],
                  Y_res,
                  Em["C3"][:, :, k],
                  R,
                  Em["n3"],
                  m3,
                  lambdax,
              )
  
          # Identify significant and insignificant terms
          Em["select"][:, k] = f_test(
              Y_id,
              Em["f0"][k],
              Y_em,
              R,
              alpha,
              m1,
              m2,
              m3,
              Em["n1"],
              Em["n2"],
              Em["n3"],
              Em["n"],
          )
  
          # Store the emulated output
          Em["Y_e"][:, k] = Em["f0"][k] + np.sum(Y_em, axis=1)
  
          # RMSE of emulator
          Em["RMSE"][k] = np.sqrt(
              np.sum(np.square(Y_id - Em["Y_e"][:, k])) / R
          )
  
          # Compute sensitivity indices
          SA["S"][:, k], SA["Sa"][:, k], SA["Sb"][:, k] = ancova(
              Y_id, Y_em, SA["V_Y"][k], R, Em["n"]
          )
  
          RT[k] = time.time() - tic  # Compute CPU time kth emulator
  
      return (SA, Em, RT, Y_em, idx)
  
imports: from typing import Dict
  from types import MethodType
  import itertools
  import time
  import warnings
  import numpy
  from numpy.linalg import solve
  from numpy.linalg import svd
  from numpy import identity
  import pandas
  from scipy import stats, special, interpolate
  from .None import common_args
  from ..util import read_param_file, ResultDict
  from SALib.plotting.hdmr import plot
  from SALib.util.problem import ProblemSpec
label_desc: Use image node_id 3 for calling a global function with example usage: _compute(X, Y, settings, init_vars) and returns: SA,  Em,  RT,  Y_em,  idx