State-space models are widely employed across various research disciplines to study unobserved dynamics. Conventional estimation techniques, such as Kalman filtering and expectation maximisation, offer valuable insights but incur high computational costs in large-scale analyses. Sparse inverse covariance estimators can mitigate these costs, but at the expense of a trade-off between enforced sparsity and increased estimation bias, necessitating careful assessment in low signal-to-noise ratio (SNR) situations. To address these challenges, we propose a three-fold solution: 1) Introducing multiple penalised state-space (MPSS) models that leverage data-driven regularisation; 2) Developing novel algorithms derived from backpropagation, gradient descent, and alternating least squares to solve MPSS models; 3) Presenting a K-fold cross-validation extension for evaluating regularisation parameters. We validate this MPSS regularisation framework through lower and more complex simulations under varying SNR conditions, including a large-scale synthetic magneto- and electro-encephalography (MEG/EEG) data analysis. In addition, we apply MPSS models to concurrently solve brain source localisation and functional connectivity problems for real event-related MEG/EEG data, encompassing thousands of sources on the cortical surface. The proposed methodology overcomes the limitations of existing approaches, such as constraints to small-scale and region-of-interest analyses. Thus, it may enable a more accurate and detailed exploration of cognitive brain functions.