Functional magnetic resonance imaging (fMRI) provides a window on the human brain at work. Spontaneous brain activity measured during resting-state has already provided many insights into brain function. In particular, recent interest in dynamic interactions between brain regions has increased the need for more advanced modeling tools. Here, we deploy a recent fMRI deconvolution technique to express resting-state temporal fluctuations as a combination of large-scale functional network activity profiles. Then, building upon a novel sparse coupled hidden Markov model (SCHMM) framework, we parameterised their temporal evolution as a mix between intrinsic dynamics, and a restricted set of cross-network modulatory couplings extracted in data-driven manner. We demonstrate and validate the method on simulated data, for which we observed that the SCHMM could accurately estimate network dynamics, revealing more precise insights about direct network-to-network modulatory influences than with conventional correlational methods. On experimental resting-state fMRI data, we unraveled a set of reproducible cross-network couplings across two independent datasets. Our framework opens new perspectives for capturing complex temporal dynamics and their changes in health and disease.