Nuclear matter is an ideal testbed for nuclear interactions with important consequences for nuclear astrophysics as well as finite nuclei. In particular, recent ab initio calculations of medium-mass to heavy nuclei have demonstrated the importance of realistic saturation properties of infinite matter for nuclear forces. We present an efficient Monte-Carlo framework for perturbative calculations of infinite nuclear matter based on two-, three-, and four-nucleon forces derived within chiral effective field theory. It enables to incorporate all many-body contributions in a transparent and also straightforward way, making it well-suited for pushing the limits of current state-of-the-art calculations to high orders in both the chiral as well as the many-body expansion. Furthermore, uncertainty estimates can be systematically extracted by order-by-order calculations, which provides important insights into the rate of convergence of each of the two expansions. After demonstrating its versatility, we make use of this novel framework to explore new chiral interactions up to next-to-next-to-next-to-leading order (N3LO) and study the equation of state of neutron and symmetric nuclear matter. Remarkably, simultaneous fits to the triton and to saturation properties can be achieved with natural 3N low-energy couplings. Taking advantage of the framework's efficacy, future chiral potentials may be optimized with respect to empirical saturation propertie