Digital Rock Physics can significantly enhance our understanding of rock behavior. However, modeling heterogeneous rocks remains challenging because of the trade-off between resolution and field of view. To address this, researchers have developed multi-scale pore network models (PNMs), which integrate PNMs from different scales to create unified multi-scale PNM. Various methodologies exist for merging PNMs from different resolutions, but they often suffer from inaccuracy, high runtime and significant memory consumption, particularly when microporosity is integrated into larger scales. This study introduces a novel fusion and an innovative upscaling approach for efficient multi-scale PNM reconstruction of rocks containing microporosity. Our methods separate resolved and unresolved porosities using different voxel sizes from CT scans at multiple resolutions. Resolved regions have larger voxel sizes, while unresolved areas retain smaller voxel sizes. We extract macroPNM from the resolved regions and generate stochastic micro-PNM for the unresolved areas. An artificial neural network (ANN), trained on micro-PNM, links micro- and macro-PNMs. The multi-scale PNMs generated using the ANN method had an average permeability of 252 +/- 3 mD, closely matching the laboratory-measured permeability of the rock (257 mD). In contrast, the average permeability of multi-scale PNMs reconstructed using the statistical method was significantly higher, at 308 +/- 38 mD. Consequently, the ANN-based reconstruction method, owing to the proper connection between scales, improved the accuracy of permeability prediction by approximately 90% compared to the statistical reconstruction method. In the next step, each microPNM is upscaled to a base pore based on its effective hydraulic conductance. These base pores are then connected to the macro-PNM using a novel approach. We utilized synchrotron CT images of an Indiana limestone rock at two resolutions as our training dataset. The single- and multi-phase flow analysis of the fused PNM demonstrated excellent agreement with laboratory-measured rock properties. Our upscaling method also reduced runtime by up to 40% (from 312 to 190 CPU-seconds) and memory consumption by approximately 68% (from 25 GB to 8 GB), all without compromising predictive accuracy.