To reduce environmental pollution and improve the efficiency of cascaded energy utilization, regional integrated energy system (RIES) has received extensive attention. An accurate multi-energy load prediction is significant for RIES as it enables stakeholders to make effective decisions for carbon peaking and carbon neutrality goals. To this end, this paper proposes a multivariate two-stage adaptive-stacking prediction (M2ASP) framework. First, a preprocessing module based on ensemble learning is proposed. The input data are preprocessed to provide a reliable database for M2ASP, and highly correlated input variables of multi-energy load prediction are determined. Then, the load prediction results of four predictors are adaptively combined in the first stage of M2ASP to enhance generalization ability. Predictor hyper-parameters and intermediate data sets of M2ASP are trained with a metaheuristic method named collaborative atomic chaotic search (CACS) to achieve the adaptive staking of M2ASP. Finally, a prediction correction of the peak load consumption period is conducted in the second stage of M2ASP. The case studies indicate that the proposed framework has higher prediction accuracy, generalization ability, and stability than other benchmark prediction models.