Dynamic models of biological systems often possess complex and multivariate mappings between input parameters and output state variables, posing challenges for comprehensive sensitivity analysis across the biologically relevant parameter space. In particular, more efficient and robust ways to obtain a solid understanding of how the sensitivity to each parameter depends on the values of the other parameters are sorely needed.We report a new methodology for global sensitivity analysis based on Hierarchical Cluster-based Partial Least Squares Regression (HC-PLSR)-based approximations (metamodelling) of the input–output mappings of dynamic models, which we expect to be generic, efficient and robust, even for systems with highly nonlinear input–output relationships. The two-step HC-PLSR metamodelling automatically separates the observations (here corresponding to different combinations of input parameter values) into groups based on the dynamic model behaviour, then analyses each group separately with Partial Least Squares Regression (PLSR). This produces one global regression model comprising all observations, as well as regional regression models within each group, where the regression coefficients can be used as sensitivity measures. Thereby a more accurate description of complex interactions between inputs to the dynamic model can be revealed through analysis of how a certain level of one input parameter affects the model sensitivity to other inputs. We illustrate the usefulness of the HC-PLSR approach on a dynamic model of a mouse heart muscle cell, and demonstrate how it reveals interaction patterns of probable biological significance not easily identifiable by a global regression-based sensitivity analysis alone.Applied for sensitivity analysis of a complex, high-dimensional dynamic model of the mouse heart muscle cell, several interactions between input parameters were identified by the two-step HC-PLSR analysis that could not be detected in the single-step global analysis. Hence, our approach has the potential to reveal new biological insight through the identification of complex parameter interaction patterns. The HC-PLSR metamodel complexity can be adjusted according to the nonlinear complexity of the input–output mapping of the analysed dynamic model through adjustment of the number of regional regression models included. This facilitates sensitivity analysis of dynamic models of varying complexities.