def predict_cdf(self, X, quantile, nsamples=200, likelihood_args=()):
        r"""
        Predictive cumulative density function of a Bayesian GLM.

        Parameters
        ----------
        X : ndarray
            (N*,d) array query input dataset (N* samples, D dimensions).
        quantile : float
            The predictive probability, :math:`p(y^* \leq \text{quantile} |
            \mathbf{x}^*, \mathbf{X}, y)`.
        nsamples : int, optional
            Number of samples for sampling the predictive CDF.
        likelihood_args : sequence, optional
            sequence of arguments to pass to the likelihood function. These are
            non-learnable parameters. They can be scalars or arrays of length
            N*.
        nsamples : int, optional
            The number of samples to draw from the posterior in order to
            approximate the predictive mean and variance.

        Returns
        -------
        p : ndarray
           The probability of y* <= quantile for the query inputs, X* of shape
           (N*,).
        p_min : ndarray
            The minimum sampled values of the predicted probability (same shape
            as p)
        p_max : ndarray
            The maximum sampled values of the predicted probability (same shape
            as p)
        """
        # Get latent function samples
        N = X.shape[0]
        ps = np.empty((N, nsamples))
        fsamples = self._sample_func(X, nsamples)

        # Push samples though likelihood cdf
        cdfarg = tuple(chain(atleast_list(self.like_hypers_), likelihood_args))
        for i, f in enumerate(fsamples):
            ps[:, i] = self.likelihood.cdf(quantile, f, *cdfarg)

        # Average transformed samples (MC integration)
        p = ps.mean(axis=1)
        p_min = ps.min(axis=1)
        p_max = ps.max(axis=1)

        return p, p_min, p_max
r"""
        Predictive cumulative density function of a Bayesian GLM.

        Parameters
        ----------
        X : ndarray
            (N*,d) array query input dataset (N* samples, D dimensions).
        quantile : float
            The predictive probability, :math:`p(y^* \leq \text{quantile} |
            \mathbf{x}^*, \mathbf{X}, y)`.
        nsamples : int, optional
            Number of samples for sampling the predictive CDF.
        likelihood_args : sequence, optional
            sequence of arguments to pass to the likelihood function. These are
            non-learnable parameters. They can be scalars or arrays of length
            N*.
        nsamples : int, optional
            The number of samples to draw from the posterior in order to
            approximate the predictive mean and variance.

        Returns
        -------
        p : ndarray
           The probability of y* <= quantile for the query inputs, X* of shape
           (N*,).
        p_min : ndarray
            The minimum sampled values of the predicted probability (same shape
            as p)
        p_max : ndarray
            The maximum sampled values of the predicted probability (same shape
            as p)