Modeling heterogeneous and fractured porous media is essential to understand the behavior of aquifers and underground reservoirs subjected to pressure and stress changes. The flow and deformation processes that occur both in the porous matrix and in the fracture network are crucial in a wide range of underground phenomena, such as land subsidence, induced seismicity, geothermal energy or green hydrogen storage. The heterogeneity of the porous matrix and the complexity of the fracture network hampers the modeling of hydromechanical phenomena in underground formations. Classical double-porosity models address the problem of upscaling coupled flow and deformation, but they can only be applied when the flux along fractures is much faster than the drainage of the porous matrix. Multirate models allow us to overcome this obstacle and address the deformation and flux processes trough fractured porous media, although their numerical implementation are still under develop. We present a new multirate formulation of multi-porosity poroelasticity that reproduces the expected behavior and scalings of coupled flow and deformation in heterogeneous fractured media. Our equivalent model describes the fractured media as a combination of several immobile regions representing the poorly impermeable matrix blocks, and a mobile zone comparable to the highly permeable fracture network. The mass exchange and the number of immobile regions considered depends on the matrix heterogeneity. To validate our formulation, we study the consolidation process of a synthetic outcrop which represents highly-heterogeneous fractured media. We compare our multirate formulation with high-fidelity numerical simulations, which explicitly model the fracture network. We use the finite element method to solve both approaches and analyze the results of land subsidence and fluid flux. The results arise that both methods agree if the characteristics of the porous matrix regions are properly chosen to represent its heterogeneity. Our multirate formulation could be used to reproduce the coupled hydromechanical behavior of the most fractured zones in reservoir-scale models, even using a coarse discretization. Therefore, subsoil fracturing could be taken into account with a low computational cost when studying fluid injection and production in deep rock formations.