def get_distribution(self, params):
        """
        Returns a scipy distribution object corresponding to the
        distribution of uncensored endog (duration) values for each
        case.

        Parameters
        ----------
        params : arrayh-like
            The model proportional hazards model parameters.

        Returns
        -------
        A list of objects of type scipy.stats.distributions.rv_discrete

        Notes
        -----
        The distributions are obtained from a simple discrete estimate
        of the survivor function that puts all mass on the observed
        failure times wihtin a stratum.
        """

        # TODO: this returns a Python list of rv_discrete objects, so
        # nothing can be vectorized.  It appears that rv_discrete does
        # not allow vectorization.

        from scipy.stats.distributions import rv_discrete

        surv = self.surv
        bhaz = self.baseline_cumulative_hazard(params)

        # The arguments to rv_discrete_float, first obtained by
        # stratum
        pk, xk = [], []

        for stx in range(self.surv.nstrat):

            exog_s = surv.exog_s[stx]

            linpred = np.dot(exog_s, params)
            if surv.offset_s is not None:
                linpred += surv.offset_s[stx]
            e_linpred = np.exp(linpred)

            # The unique failure times for this stratum (the support
            # of the distribution).
            pts = bhaz[stx][0]

            # The individual cumulative hazards for everyone in this
            # stratum.
            ichaz = np.outer(e_linpred, bhaz[stx][1])

            # The individual survival functions.
            usurv = np.exp(-ichaz)
            usurv = np.concatenate((usurv, np.zeros((usurv.shape[0], 1))),
                                   axis=1)

            # The individual survival probability masses.
            probs = -np.diff(usurv, 1)

            pk.append(probs)
            xk.append(np.outer(np.ones(probs.shape[0]), pts))

        # Pad to make all strata have the same shape
        mxc = max([x.shape[1] for x in xk])
        for k in range(self.surv.nstrat):
            if xk[k].shape[1] < mxc:
                xk1 = np.zeros((xk.shape[0], mxc))
                pk1 = np.zeros((pk.shape[0], mxc))
                xk1[:, -mxc:] = xk
                pk1[:, -mxc:] = pk
                xk[k], pk[k] = xk1, pk1

        xka = np.nan * np.zeros((len(self.endog), mxc), dtype=np.float64)
        pka = mxc / np.ones((len(self.endog), mxc), dtype=np.float64)

        for stx in range(self.surv.nstrat):
            ix = self.surv.stratum_rows[stx]
            xka[ix, :] = xk[stx]
            pka[ix, :] = pk[stx]

        dist = rv_discrete_float(xka, pka)

        return dist
dataset/ETHPy150Open statsmodels/statsmodels/statsmodels/duration/hazard_regression.py/PHReg.get_distribution/swapped operands of `/`