This paper formulates and addresses the problem of equivalence in terms of multistability properties between nonlinear models of gene regulatory systems of different dimensionality. Given a nonlinear dynamical model of a gene regulatory network and the structure of another higher-dimensional gene regulatory network, the aim is to find a dynamical model for the latter that has the same equilibria and stability properties as the former. We propose construction rules for the dynamics of a high-dimensional system, given the low-dimensional system and the high-dimensional network structure. These construction rules yield a multistability-equivalent system, as we prove in this work. We demonstrate the value of our method by applying it to an example of a multistable gene regulatory network involved in mesenchymal stem cell differentiation. Here, differentiation is described by a core motif of three genetic regulators, but the detailed network contains at least nine genes. The proposed construction method allows to transfer the multistability based differentation mechanism of the core motif to the more detailed gene regulatory network. Copyright MULTISTABILITY EQUIVALENCE BETWEEN GENE REGULATORY NETWORKS 4149 A central property of GRNs is that many of them exhibit multistability, thereby realizing the coexistence of several stable states of the system, of which a suitable one can be adopted under the respective biological circumstances ([2, 5, 7, 8]). The analysis of multistability properties of such systems is of high interest and has attracted the development of new specially suited methods [9][10][11][12].It is desired that the analysis of a conceptual, low-dimensional model should allow for conclusions about a higher-dimensional model. On the one hand, low-dimensional systems are amenable to systems-theoretic tools such as multistability and bifurcation analysis. On the other hand, a high-dimensional dynamic model might be required for fitting model dynamics to data from gene expression measurements, and thus for quantitative predictions about the real biological system. From this, the task arises how a dynamic model for a high-dimensional system can be obtained if only a core motif model of the gene regulatory dynamics, implemented as a low-dimensional ODE system, is available. Often, the structure of the high-dimensional network is known from network identification methods (or at least a small number of alternatives exist), but the dynamics in form of an ODE system are unknown and need to be constructed. The high-dimensional ODE system to be obtained should have the same number of equilibria, and each equilibrium should have the same stability properties, as the low-dimensional system. This is somewhat complementary to the fields of model reduction and model decomposition: In model reduction, one aims to reduce a high-dimensional model to a low-dimensional one [13]. The aim of model decomposition, in contrast, is to partition a high-dimensional model into weakly connected modules, which can be analyzed by systems-...