As all kinds of physics-based and data-driven models are emerging in hydrologic and hydraulic engineering, Bayesian model averaging (BMA) is one of the popular multi-model methods used to account for various uncertainty sources in the flood modeling process and generate robust ensemble predictions. The reliability of BMA parameters (weights and variances) determines the accuracy of BMA predictions. However, the uncertainty in BMA parameters with fixed values, which are usually obtained from Expectation-Maximization (EM) algorithm, has not been adequately investigated in BMA-related applications over the past few decades. Given the limitations of the commonly used EM algorithm, Metropolis-Hastings (M-H) algorithm, which is one of the most widely used algorithms in Markov Chain Monte Carlo (MCMC) method, is proposed to estimate BMA parameters. Both numerical experiments and one-dimensional HEC-RAS models are employed to examine the applicability of M-H algorithm with multiple independent Markov chains. The performances of EM and M-H algorithms are compared based on the daily water stage predictions from 10 model members. Results show that BMA weights estimated from both algorithms are comparable, while BMA variances obtained from M-H algorithm are closer to the given variances in the numerical experiment. Moreover, the normal proposal used in M-H algorithm can yield narrower distributions for BMA weights than those from the uniform proposal. Overall, MCMC approach with multiple chains can provide more information associated with the uncertainty of BMA parameters and its performance is better than the default EM algorithm in terms of multiple evaluation metrics as well as algorithm flexibility.