Distribution and connectivity of pores in shale is difficult to characterise due to their nano-micro multiscale characteristics. Flow in the porous structure is correspondingly controlled by multiple mechanisms and thus it is challenging to estimate the flow properties of shale gas. In view of this, the pore network was employed in this study to represent the porous structure and carry out flow simulation in an organic-rich shale matrix. A self-similarity based stochastic method, which combines the advantages of μ-CT and SEM, was used to characterise the porous structure with a resolution of 10 nm. A pore network of size of 10 × 10 × 10 μm3 was constructed by the maximum ball algorithm and the multi-physics flow (including the wide range non-continuum bulk flow and surface diffusion) was simulated in the reconstructed network. The simulation results show that the apparent gas permeability increases by a factor of 3.16, with 36% contribution arising from surface diffusion when the downstream pressure depletes from 9 MPa to 2 MPa. Sensitivity analyses imply that the apparent gas permeability is dependent on the size and shape of throats, compressibility factor of gas, Langmuir adsorption parameters and reservoir conditions.