Focusing on reliability analysis of nonparametric probability boxes (p-boxes), an adaptive discrete cumulative distribution function method based on Bayesian active learning is proposed to efficiently estimate failure probability bounds. First, a Bayesian active learning framework is employed to train a Gaussian Process surrogate model, replacing direct functional evaluations in the probability integration process. Second, an explicit expression for the error induced by discretization grid spacing is derived, establishing a quantitative estimation of the numerical integration error bound. Third, an adaptive discretization criterion based on a comprehensive error estimator is developed; this criterion allocates new nodes to subintervals that most significantly impact the accuracy of failure probability bounds, achieving a locally optimal allocation of computational resources. Finally, under a convergence criterion controlled by the upper bound of integration error, convergent upper and lower bounds of the failure probability are obtained. Numerical examples demonstrate the effectiveness and efficiency of the proposed method.