A probabilistic load flow (PLF) is an effective tool that helps describe the uncertainty of power system operation. However, when confronting random variables with non-Gaussian distributions and highly discrete characteristics, existing PLF methods have difficulty balancing efficiency and accuracy. Therefore, a novel approach based on bivariate dimension reduction (BDR) and the Johnson system is proposed herein. BDR is used to estimate the moments of output random variables (ORVs). Because BDR considers the joint effect of input random variables, it significantly reduces the estimation error for high-order moments in particular. In addition, a strategy to improve BDR efficiency is proposed. The Johnson system is used to obtain the probability distributions of ORVs as it has better adaptability and accuracy than the series expansion method. Case studies including comparisons between this approach and others found in the literature were conducted, and the results obtained showed that the proposed method has better performance than previous approaches. INDEX TERMS Bivariate dimension reduction, Johnson system, probabilistic load flow, renewable energy source.