code: def mediation_analysis(
      data=None,
      x=None,
      m=None,
      y=None,
      covar=None,
      alpha=0.05,
      n_boot=500,
      seed=None,
      return_dist=False,
      logreg_kwargs=None,
  ):
      """Mediation analysis using a bias-correct non-parametric bootstrap method.
  
      Parameters
      ----------
      data : :py:class:`pandas.DataFrame`
          Dataframe.
      x : str
          Column name in data containing the predictor variable.
          The predictor variable must be continuous.
      m : str or list of str
          Column name(s) in data containing the mediator variable(s).
          The mediator(s) can be continuous or binary (e.g. 0 or 1).
          This function supports multiple parallel mediators.
      y : str
          Column name in data containing the outcome variable.
          The outcome variable must be continuous.
      covar : None, str, or list
          Covariate(s). If not None, the specified covariate(s) will be included
          in all regressions.
      alpha : float
          Significance threshold. Used to determine the confidence interval,
          :math:`\\text{CI} = [\\alpha / 2 ; 1 - \\alpha / 2]`.
      n_boot : int
          Number of bootstrap iterations for confidence intervals and p-values
          estimation. The greater, the slower.
      seed : int or None
          Random state seed.
      logreg_kwargs : dict or None
          Dictionary with optional arguments passed to :py:func:`pingouin.logistic_regression`
      return_dist : bool
          If True, the function also returns the indirect bootstrapped beta
          samples (size = n_boot). Can be plotted for instance using
          :py:func:`seaborn.distplot()` or :py:func:`seaborn.kdeplot()`
          functions.
  
      Returns
      -------
      stats : :py:class:`pandas.DataFrame`
          Mediation summary:
  
          * ``'path'``: regression model
          * ``'coef'``: regression estimates
          * ``'se'``: standard error
          * ``'CI[2.5%]'``: lower confidence interval
          * ``'CI[97.5%]'``: upper confidence interval
          * ``'pval'``: two-sided p-values
          * ``'sig'``: statistical significance
  
      See also
      --------
      linear_regression, logistic_regression
  
      Notes
      -----
      Mediation analysis [1]_ is a *"statistical procedure to test
      whether the effect of an independent variable X on a dependent variable
      Y (i.e., X → Y) is at least partly explained by a chain of effects of the
      independent variable on an intervening mediator variable M and of the
      intervening variable on the dependent variable (i.e., X → M → Y)"* [2]_.
  
      The **indirect effect** (also referred to as average causal mediation
      effect or ACME) of X on Y through mediator M quantifies the estimated
      difference in Y resulting from a one-unit change in X through a sequence of
      causal steps in which X affects M, which in turn affects Y.
      It is considered significant if the specified confidence interval does not
      include 0. The path 'X --> Y' is the sum of both the indirect and direct
      effect. It is sometimes referred to as total effect.
  
      A linear regression is used if the mediator variable is continuous and a
      logistic regression if the mediator variable is dichotomous (binary).
      Multiple parallel mediators are also supported.
  
      This function will only work well if the outcome variable is continuous.
      It does not support binary or ordinal outcome variable. For more
      advanced mediation models, please refer to the
      `lavaan <http://lavaan.ugent.be/tutorial/mediation.html>`_ or  `mediation
      <https://cran.r-project.org/web/packages/mediation/mediation.pdf>`_ R
      packages, or the `PROCESS macro
      <https://www.processmacro.org/index.html>`_ for SPSS.
  
      The two-sided p-value of the indirect effect is computed using the
      bootstrap distribution, as in the mediation R package. However, the p-value
      should be interpreted with caution since it is not constructed
      conditioned on a true null hypothesis [3]_ and varies depending on the
      number of bootstrap samples and the random seed.
  
      Note that rows with missing values are automatically removed.
  
      Results have been tested against the R mediation package and this tutorial
      https://data.library.virginia.edu/introduction-to-mediation-analysis/
  
      References
      ----------
      .. [1] Baron, R. M. & Kenny, D. A. The moderator–mediator variable
             distinction in social psychological research: Conceptual, strategic,
             and statistical considerations. J. Pers. Soc. Psychol. 51, 1173–1182
             (1986).
  
      .. [2] Fiedler, K., Schott, M. & Meiser, T. What mediation analysis can
             (not) do. J. Exp. Soc. Psychol. 47, 1231–1236 (2011).
  
  
      .. [3] Hayes, A. F. & Rockwood, N. J. Regression-based statistical
             mediation and moderation analysis in clinical research:
             Observations, recommendations, and implementation. Behav. Res.
             Ther. 98, 39–57 (2017).
  
      Code originally adapted from https://github.com/rmill040/pymediation.
  
      Examples
      --------
      1. Simple mediation analysis
  
      >>> from pingouin import mediation_analysis, read_dataset
      >>> df = read_dataset('mediation')
      >>> mediation_analysis(data=df, x='X', m='M', y='Y', alpha=0.05,
      ...                    seed=42)
             path      coef        se          pval  CI[2.5%]  CI[97.5%]  sig
      0     M ~ X  0.561015  0.094480  4.391362e-08  0.373522   0.748509  Yes
      1     Y ~ M  0.654173  0.085831  1.612674e-11  0.483844   0.824501  Yes
      2     Total  0.396126  0.111160  5.671128e-04  0.175533   0.616719  Yes
      3    Direct  0.039604  0.109648  7.187429e-01 -0.178018   0.257226   No
      4  Indirect  0.356522  0.083313  0.000000e+00  0.219818   0.537654  Yes
  
      2. Return the indirect bootstrapped beta coefficients
  
      >>> stats, dist = mediation_analysis(data=df, x='X', m='M', y='Y',
      ...                                  return_dist=True)
      >>> print(dist.shape)
      (500,)
  
      3. Mediation analysis with a binary mediator variable
  
      >>> mediation_analysis(data=df, x='X', m='Mbin', y='Y', seed=42).round(3)
             path   coef     se   pval  CI[2.5%]  CI[97.5%]  sig
      0  Mbin ~ X -0.021  0.116  0.857    -0.248      0.206   No
      1  Y ~ Mbin -0.135  0.412  0.743    -0.952      0.682   No
      2     Total  0.396  0.111  0.001     0.176      0.617  Yes
      3    Direct  0.396  0.112  0.001     0.174      0.617  Yes
      4  Indirect  0.002  0.050  0.960    -0.072      0.146   No
  
      4. Mediation analysis with covariates
  
      >>> mediation_analysis(data=df, x='X', m='M', y='Y',
      ...                    covar=['Mbin', 'Ybin'], seed=42).round(3)
             path   coef     se   pval  CI[2.5%]  CI[97.5%]  sig
      0     M ~ X  0.559  0.097  0.000     0.367      0.752  Yes
      1     Y ~ M  0.666  0.086  0.000     0.495      0.837  Yes
      2     Total  0.420  0.113  0.000     0.196      0.645  Yes
      3    Direct  0.064  0.110  0.561    -0.155      0.284   No
      4  Indirect  0.356  0.086  0.000     0.209      0.553  Yes
  
      5. Mediation analysis with multiple parallel mediators
  
      >>> mediation_analysis(data=df, x='X', m=['M', 'Mbin'], y='Y',
      ...                    seed=42).round(3)
                  path   coef     se   pval  CI[2.5%]  CI[97.5%]  sig
      0          M ~ X  0.561  0.094  0.000     0.374      0.749  Yes
      1       Mbin ~ X -0.005  0.029  0.859    -0.063      0.052   No
      2          Y ~ M  0.654  0.086  0.000     0.482      0.825  Yes
      3       Y ~ Mbin -0.064  0.328  0.846    -0.715      0.587   No
      4          Total  0.396  0.111  0.001     0.176      0.617  Yes
      5         Direct  0.040  0.110  0.721    -0.179      0.258   No
      6     Indirect M  0.356  0.085  0.000     0.215      0.538  Yes
      7  Indirect Mbin  0.000  0.010  0.952    -0.017      0.025   No
      """
      # Sanity check
      assert isinstance(x, (str, int)), "y must be a string or int."
      assert isinstance(y, (str, int)), "y must be a string or int."
      assert isinstance(
          m, (list, str, int)
      ), "Mediator(s) must be a list, string or int."
      assert isinstance(covar, (type(None), str, list, int))
      if isinstance(m, (str, int)):
          m = [m]
      n_mediator = len(m)
      assert isinstance(data, pd.DataFrame), "Data must be a DataFrame."
      # Check for duplicates
      assert n_mediator == len(
          set(m)
      ), "Cannot have duplicates mediators."
      if isinstance(covar, (str, int)):
          covar = [covar]
      if isinstance(covar, list):
          assert len(covar) == len(
              set(covar)
          ), "Cannot have duplicates covar."
          assert set(m).isdisjoint(
              covar
          ), "Mediator cannot be in covar."
      # Check that columns are in dataframe
      columns = _fl([x, m, y, covar])
      keys = data.columns
      assert all(
          [c in keys for c in columns]
      ), "Column(s) are not in DataFrame."
      # Check that columns are numeric
      err_msg = "Columns must be numeric or boolean."
      assert all(
          [data[c].dtype.kind in "bfiu" for c in columns]
      ), err_msg
  
      # Drop rows with NAN Values
      data = data[columns].dropna()
      n = data.shape[0]
      assert n > 5, "DataFrame must have at least 5 samples (rows)."
  
      # Check if mediator is binary
      mtype = "logistic" if all(data[m].nunique() == 2) else "linear"
  
      # Check if a dict with kwargs for logistic_regression has been passed
      logreg_kwargs = {} if logreg_kwargs is None else logreg_kwargs
  
      # Name of CI
      ll_name = "CI[%.1f%%]" % (100 * alpha / 2)
      ul_name = "CI[%.1f%%]" % (100 * (1 - alpha / 2))
  
      # Compute regressions
      cols = ["names", "coef", "se", "pval", ll_name, ul_name]
  
      # For speed, we pass np.array instead of pandas DataFrame
      X_val = data[
          _fl([x, covar])
      ].to_numpy()  # X + covar as predictors
      XM_val = data[
          _fl([x, m, covar])
      ].to_numpy()  # X + M + covar as predictors
      M_val = data[m].to_numpy()  # M as target (no covariates)
      y_val = data[y].to_numpy()  # y as target (no covariates)
  
      # For max precision, make sure rounding is disabled
      old_options = options.copy()
      options["round"] = None
  
      # M(j) ~ X + covar
      sxm = {}
      for idx, j in enumerate(m):
          if mtype == "linear":
              sxm[j] = linear_regression(
                  X_val, M_val[:, idx], alpha=alpha
              ).loc[[1], cols]
          else:
              sxm[j] = logistic_regression(
                  X_val, M_val[:, idx], alpha=alpha, **logreg_kwargs
              ).loc[[1], cols]
          sxm[j].at[1, "names"] = "%s ~ X" % j
      sxm = pd.concat(sxm, ignore_index=True)
  
      # Y ~ M + covar
      smy = linear_regression(
          data[_fl([m, covar])], y_val, alpha=alpha
      ).loc[1:n_mediator, cols]
      # Average Total Effects (Y ~ X + covar)
      sxy = linear_regression(X_val, y_val, alpha=alpha).loc[[1], cols]
      # Average Direct Effects (Y ~ X + M + covar)
      direct = linear_regression(XM_val, y_val, alpha=alpha).loc[
          [1], cols
      ]
  
      # Rename paths
      smy["names"] = smy["names"].apply(lambda x: "Y ~ %s" % x)
      direct.at[1, "names"] = "Direct"
      sxy.at[1, "names"] = "Total"
  
      # Concatenate and create sig column
      stats = pd.concat((sxm, smy, sxy, direct), ignore_index=True)
      stats["sig"] = np.where(stats["pval"] < alpha, "Yes", "No")
  
      # Bootstrap confidence intervals
      rng = np.random.RandomState(seed)
      idx = rng.choice(np.arange(n), replace=True, size=(n_boot, n))
      ab_estimates = np.zeros(shape=(n_boot, n_mediator))
      for i in range(n_boot):
          ab_estimates[i, :] = _point_estimate(
              X_val,
              XM_val,
              M_val,
              y_val,
              idx[i, :],
              n_mediator,
              mtype,
              **logreg_kwargs,
          )
  
      ab = _point_estimate(
          X_val,
          XM_val,
          M_val,
          y_val,
          np.arange(n),
          n_mediator,
          mtype,
          **logreg_kwargs,
      )
      indirect = {
          "names": m,
          "coef": ab,
          "se": ab_estimates.std(ddof=1, axis=0),
          "pval": [],
          ll_name: [],
          ul_name: [],
          "sig": [],
      }
  
      for j in range(n_mediator):
          ci_j = _bias_corrected_ci(
              ab_estimates[:, j], indirect["coef"][j], alpha=alpha
          )
          indirect[ll_name].append(min(ci_j))
          indirect[ul_name].append(max(ci_j))
          # Bootstrapped p-value of indirect effect
          # Note that this is less accurate than a permutation test because the
          # bootstrap distribution is not conditioned on a true null hypothesis.
          # For more details see Hayes and Rockwood 2017
          indirect["pval"].append(
              _pval_from_bootci(ab_estimates[:, j], indirect["coef"][j])
          )
          indirect["sig"].append(
              "Yes" if indirect["pval"][j] < alpha else "No"
          )
  
      # Create output dataframe
      indirect = pd.DataFrame.from_dict(indirect)
      if n_mediator == 1:
          indirect["names"] = "Indirect"
      else:
          indirect["names"] = indirect["names"].apply(
              lambda x: "Indirect %s" % x
          )
      stats = pd.concat(
          [stats, indirect], axis=0, ignore_index=True, sort=False
      )
      stats = stats.rename(columns={"names": "path"})
  
      # Restore options
      options.update(old_options)
  
      if return_dist:
          return _postprocess_dataframe(stats), np.squeeze(ab_estimates)
      else:
          return _postprocess_dataframe(stats)
  
file_path: pingouin/pingouin/regression.py
imports: import itertools
  import warnings
  import numpy
  import pandas
  import pandas_flavor
  from scipy.stats import t, norm
  from scipy.linalg import pinvh, lstsq
  from pingouin.config import options
  from pingouin.utils import remove_na
  from pingouin.utils import _flatten_list
  from pingouin.utils import _postprocess_dataframe
label_desc: Use image node_id 7 for calling a global function with example usage: mediation_analysis(data, x, m, y, covar, alpha, n_boot, seed, return_dist, logreg_kwargs) and returns: _postprocess_dataframe,  np, _postprocess_dataframe