def fit(self, method="parametric", n_rep=1000):
        """
        Fit a regression model to assess mediation.

        Arguments
        ---------
        method : string
            Either 'parametric' or 'bootstrap'.
        n_rep : integer
            The number of simulation replications.

        Returns a MediationResults object.
        """

        if method.startswith("para"):
            # Initial fit to unperturbed data.
            outcome_result = self._fit_model(self.outcome_model, self._outcome_fit_kwargs)
            mediator_result = self._fit_model(self.mediator_model, self._mediator_fit_kwargs)
        elif not method.startswith("boot"):
            raise("method must be either 'parametric' or 'bootstrap'")

        indirect_effects = [[], []]
        direct_effects = [[], []]

        for iter in range(n_rep):

            if method == "parametric":
                # Realization of outcome model parameters from sampling distribution
                outcome_params = self._simulate_params(outcome_result)

                # Realization of mediation model parameters from sampling distribution
                mediation_params = self._simulate_params(mediator_result)
            else:
                outcome_result = self._fit_model(self.outcome_model,
                                                 self._outcome_fit_kwargs, boot=True)
                outcome_params = outcome_result.params
                mediator_result = self._fit_model(self.mediator_model,
                                                  self._mediator_fit_kwargs, boot=True)
                mediation_params = mediator_result.params

            # predicted outcomes[tm][te] is the outcome when the
            # mediator is set to tm and the outcome/exposure is set to
            # te.
            predicted_outcomes = [[None, None], [None, None]]
            for tm in 0, 1:
                mex = self._get_mediator_exog(tm)
                gen = self.mediator_model.get_distribution(mediation_params,
                                                           mediator_result.scale,
                                                           exog=mex)
                potential_mediator = gen.rvs(mex.shape[0])

                for te in 0, 1:
                    oex = self._get_outcome_exog(te, potential_mediator)
                    po = self.outcome_model.predict(outcome_params, oex)
                    predicted_outcomes[tm][te] = po

            for t in 0, 1:
                indirect_effects[t].append(predicted_outcomes[1][t] - predicted_outcomes[0][t])
                direct_effects[t].append(predicted_outcomes[t][0] - predicted_outcomes[t][1])

        for t in 0, 1:
            indirect_effects[t] = np.asarray(indirect_effects[t]).T
            direct_effects[t] = np.asarray(direct_effects[t]).T

        self.indirect_effects = indirect_effects
        self.direct_effects = direct_effects

        rslt = MediationResults(self.indirect_effects, self.direct_effects)
        rslt.method = method
        return rslt
dataset/ETHPy150Open statsmodels/statsmodels/statsmodels/stats/mediation.py/Mediation.fit/swapped operands of `-`