  The regulation of the cell state is a complex process involving several
components. These complex dynamics can be modeled using Boolean networks,
allowing us to explain the existence of different cell states and the
transition between them. Boolean models have been introduced both as specific
examples and as ensemble or distribution network models. However, current
ensemble Boolean network models do not make a systematic distinction between
different cell components such as epigenetic factors, gene and transcription
factors. Consequently, we still do not understand their relative contributions
in controlling the cell fate. In this work we introduce and study higher order
Boolean networks, which feature an explicit distinction between the different
cell components and the types of interactions between them. We show that the
stability of the cell state dynamics can be determined solving the eigenvalue
problem of a matrix representing the regulatory interactions and their
strengths. The qualitative analysis of this problem indicates that, in addition
to the classification into stable and chaotic regimes, the cell state can be
simple or complex depending on whether it can be deduced from the independent
study of its components or not. Finally, we illustrate how the model can be
expanded considering higher levels and higher order dynamics.
