{'fraction_non_alphanumeric': 0.023956984964651998, 'fraction_numerical': 0.006770885193667231, 'mean_word_length': 4.944832484905883, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 10, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Higher-order percolation processes on multiplex hypergraphs

Higher order interactions are increasingly recognised as a fundamental aspect of complex systems ranging from the brain to social contact networks. Hypergraph as well as simplicial complexes capture the higher-order interactions of complex systems and allow to investigate the relation between their higher-order structure and their function. Here we establish a general framework for assessing hypergraph robustness and we characterize the critical properties of simple and higher-order percolation processes. This general framework builds on the formulation of the random multiplex hypergraph ensemble where each layer is characterized by hyperedges of given cardinality. We reveal the relation between higher-order percolation processes in random multiplex hypergraphs, interdependent percolation of multiplex networks and K-core percolation. The structural correlations of the random multiplex hypergraphs are shown to have a significant effect on their percolation properties. The wide range of critical behaviors observed for higher-order percolation processes on multiplex hypergraphs elucidates the mechanisms responsible for the emergence of discontinuous transition and uncovers interesting critical properties which can be applied to the study of epidemic spreading and contagion processes on higher-order networks.

Higher order interactions are increasingly recognised as a fundamental aspect of complex systems ranging from the brain to social contact networks. Hypergraph as well as simplicial complexes capture the higher-order interactions of complex systems and allow to investigate the relation between their higher-order structure and their function. Here we establish a general framework for assessing hypergraph robustness and we characterize the critical properties of simple and higher-order percolation processes. This general framework builds on the formulation of the random multiplex hypergraph ensemble where each layer is characterized by hyperedges of given cardinality. We reveal the relation between higher-order percolation processes in random multiplex hypergraphs, interdependent percolation of multiplex networks and K-core percolation. The structural correlations of the random multiplex hypergraphs are shown to have a significant effect on their percolation properties. The wide range of critical behaviors observed for higher-order percolation processes on multiplex hypergraphs elucidates the mechanisms responsible for the emergence of discontinuous transition and uncovers interesting critical properties which can be applied to the study of epidemic spreading and contagion processes on higherorder networks.
Higher-order networks include both hypergraphs and simplicial complexes and encode the set of higher-order interactions present in systems as different as social [9][10][11][12], ecological [13] and brain networks [14]. Multilayer networks represent complex systems in which interactions of different nature and connotation can exist forming networks of networks. As such multilayer networks and in particular multiplex networks are becoming the new paradigm to describe social, financial as well as biological networks [6][7][8].
In this work we investigate the interplay between structure and the dynamics of higher-order networks providing a comprehensive multilayer framework to study higher-order percolation processes on hypergraphs.
Percolation [29][30][31][32] is a fundamental dynamical process defined on networks that predicts the fraction of nodes in the giant component of a network. Having a non-zero giant component is the minimal requisite for observing collective phenomena on networks, emerging from epidemic spreading, diffusion and opinion dynamics. Therefore studying percolation of a given network has important consequences for investigating a wide range of dynamical properties defined on networks.
Percolation theory has been extensively studied in single networks since the early days of Network Science [29][30][31][32]. In particular node and link percolation have been investigated in random networks with arbitrary degree distribution. In node percolation nodes are initially damaged with probability 1− p, in link percolation links are initially damaged with probability 1 − p. In both percolation models the fraction of nodes in the giant component is studied as a function of p characterizing how the network is dismantled/disconnected by increasing the entity of the initial damage. Interestingly it has been found that the critical properties of percolation are strongly affected by the network topology of the underlying network. In particular a classic result of percolation theory is that scale-free networks are robust to random damage and in the infinite network limit can sustain a non-zero fraction of nodes in the giant component for any non-zero value of p [33,34]. While percolation on single networks is a continuous second order transition, K-core percolation [35,36], studying the emergence of the K-core with K ≥ 2 in complex networks can display a discontinuous hybrid transition if the degree distribution of the network has finite second moment.
With the recent surge of interest on generalized network structures percolation theory has further expanded thanks to the formulation of the interdependent percolation in multiplex networks that displays a discontinuous hybrid transition in correspondence of large avalanches of failure events [6,[37][38][39]. In interdependent percolation the order parameter is the fraction of nodes in the mutually connected giant component that is the giant component formed by nodes connected by at least one path in each layer of the multiplex network. Interdependent percolation on multiplex networks is highly affected by the correlations [15,16,40] of the underlying multiplex network structure. Indeed both interlayer degree correlations [15] and link overlap [41][42][43] have been shown to have a very significant effect on the critical properties of interdependent percolation. This field has been growing at a very fast pace and many results related to the robustness and resilience of multiplex networks have been obtained including the formulation of interdependent percolation in network of networks [44], weak percolation [45,46], optimal percolation [47,48], combinatorial optimization problems [49], K-core multiplex percolation [50] and percolation with redundant interdependencies [51].
This very important subject in network theory has contributed to a much deeper theoretical understanding of the mechanisms leading to discontinuous percolation transitions in complex networks (see recent review articles [52,53]).
Recently the rich interplay between network geometry and topology and the critical properties of the percolation transition on higher-order networks has also gained increasing attention. Hyperbolic simplicial complexes which can be treated within the real-space renormalization group [54][55][56][57][58], have been shown to reveal a rich phase diagram including discontinuous transitions for standard link percolation. Moreover homological percolation [59,60] has been show to characterize the emergence of a non trivial homology for higher-order network topologies.
If we focus on random hypergraphs which do not display an hyperbolic network geometry, the critical properties of percolation an higher-order percolation process are just starting to be explored. Indeed some interesting results related to core percolation in hypergraphs have been recently published in Ref. [61]. However so far the investigation of percolation on random hypegraphs has been restricted to very simple cases of hypergraphs which hyperedges have fixed cardinality [62,63].
In this paper we relate higher-order percolation on hypergraphs to generalized percolation processes in multiplex networks. Random hypergraphs can have a non-trivial underlying multiplex topology in which each layer capture the set of hyperedges of a given cardinality. This allows us to define ensembles of random multiplex hypergraphs in which each node is assigned set of generalized degrees {k [m] i } where k [m] i indicates the number of hyperedges of cardinality m incident to node i. As such multiplex hypergraphs are characterized by important interlayer generalized hyperdegree correlations.
Here we show that standard percolation is affected by the nontrivial topology of multiplex hypergraphs and by their interlayer correlations that can be tuned to increase of decrease the percolation threshold of the hypergraph. Most importantly our work reveal how the multiplex nature of the multiplex hypergraph ensembles can be exploited to propose higher-order percolation problems displaying a rich interplay between higherorder topology and dynamics and a rich set of phenomena, including discontinuous hybrid transitions and multiple percolation transitions.
The paper is structured as follows: in Sec II we present the random multiplex hypergraph model and we compare the model with the already widely used model of random hypergraphs; in Sec III and Sec. IV we investigate the properties of standard node and link percolation on the random hypergraphs and on the random multiplex hypegraphs respectively; in Sec. V we provide a general framework to study higher-order percolation processes on random multiplex hypergraphs; finally in Sec. VI we provide the concluding remarks.

A. Random hypergraphs
In this paragraph we introduce random hypergraphs used widely in the literature. This model will be subsequently compared with the model of random multiplex hypergraph which allows us to capture more rich hypegraphs topologies. Hypegraphs H = (V, H) are formed by a set V of N nodes and a set H of hyperedges of different cardinality m ≤ M. The number of hyperedges incident to a node is also called its hyperdegree. Therefore if all hyperedges have cardinality m = 2, i.e. all hyperedges are essentially links describing pairwise interactions, then the hypergraph reduce to a network, and the definition of hyperdegree reduces to the definition of degree. In general, in hypergraphs containing hyperedges of different cardinality, the hyperdegree counts the number of hyperedges incident to a node regardless of their cardinality.
The simplest model of hypergraph here called the random hypergraph model is a maximum entropy hypergraph model with given hyperdegree distribution P(k) and distribu-tionP(m) of hyperedge cardinalities. Therefore as far as a node has a given hyperdegree k drawn from the hyperdegree distribution P(k), the model is agnostic on the cardinality of its incident hyperedges.
The maximum entropy ensemble of random hypergraphs is a simple model that is strictly related to maximum entropy factor graphs. Factor graphs are bipartite networks G B (V, U, E) formed by a set of nodes V and a set of factor nodes U which do not overlap and a set E of pairwise interactions with each interaction linking a node to a factor node. The factor graph is simply related to the hypergraph by a simple mapping. The set V of nodes of the factor graph maps to the set V nodes of the hypergraph. Each factor node of the set U is in one-to-one correspondence with the hyperedges of the hypergraph. The hyperdegree distribution P(k) and distribu-tionP(m) of hyperedge cardinalities correspond to the degree distribution of the nodes and of the factor nodes of the factor graph respectively.
In the uncorrelated hypergraph ensemble, as long as the hyperdegree distribution and the distribution of hyperegde cardinality have a structural cutoff, the probability that node i is connnected to hyperedge α is given bỹ where k i indicates the hyperdegree of node i and m α indicates the cardinality of the hyperedge/factor node α. The corresponding hypergraph includes and hyperedge α =

B. Random multiplex hypergraphs
The random hypergraphs described in the previous paragraph are maximum entropy ensemble in which we fix the hyperdegree of each node. Here we consider hypergraphs in which we assign to each node a set of generalized hyperdegrees each one fixing the incident number of hyperedges of a given cardinality. This allow to control the number of hyperedges of a given cardinality incident to each node and providing a more refined hypergraph model than the random hypergraph.
As we will see this model can be mapped to a multiplex network model [6,40] hence we indicate this model as random multiplex hypergraph.
In this case, as in the previous case, we consider an ensem- where k [m] i indicates the number of hyperedges of degree m incident node i, i.e. , For these hypergraphs we can define the generalized hyperdegree distribution P(k) as the probability that a random node of the hypergraph has generalized hyperdegrees The random multiplex hypergraph is the maximum entropy hypergraph model with given generalized hyperdegree distribution P(k) and given distributionP(m) of cardinality of the hyperedges.
In the hypergraph setting there are no constraints relating hyperedges of different cardinality. Therefore the hyperedges of different cardinality can be drawn independently. This is different from what happens in the simplicial complex setting where if a simplex indicating interaction between m nodes is present in the simplicial complex, also all the simplices formed by the proper subset of its nodes belong to the simplicial complex. Despite this difference, simplicial complex models can be very efficiently used to model hypergraphs. Indeed the maximum entropy hypergraph with given generalized hyperdegree sequences can be constructed starting from the well establish configuration model of pure simplicial complexes [4] by mapping the hypergraph to a multiplex network [12] in which every layer indicates the interactions described by hyperedges of a given size m.
The algorithm to construct a random multiplex hypergraph (see Figure 1) is: (1) Consider a multiplex networks with M −1 layers m with 2 ≤ m < M and N nodes corresponding to the N nodes of the hypergraph.
(2) For each layer m consider configuration model of pure (m − 1)-dimensional simplicial complexes [4] with generalized degree sequence In this hypothesis the probability p [m] [i 1 ,i 2 ...,i m ] of the hyperedge [i 1 , i 2 . . . , i m ] is given by [4] Note that some layers might be empty if In this case the number of layers of the random multiplex hypegraph is given by the number M of layers with at least one hyperedge.
(3) Consider the hypergraph obtained by aggregating all the layers, i.e. considering all the interactions of different sizes 2 ≤ m < M. Note that this aggregated hypergraph, differently from the aggregated multiplex network with pairwise interactions, retains its multilayer nature as the hyperedges of different cadinality can be easily distinguished also in the aggregated version of the hypergraph. Therefore we will not make a distinction between this aggregated hypergraphs and their multiplex representation.
The random multiplex hypergraph ensemble can also be viewed as a multiplex networks of factor graphs where each layer is a factor graph in which factor nodes have fixed degree. In this factor graph interpretation, the probability that node i is connected to factor node α in a given layer describing m-body interactions, i.e. with m α = m is given bỹ The mapping between the hypergraph model with given generalized hyperdegree sequences and mutliplex networks, allows us to address the role of correlations have in this hypergraph model. Indeed, by considering a parallelism to multiplex networks we can investigate different types of possible correlations in hypergraph models. First of all the hyperdegrees k of a given nodes can be correlated. In a hypergraph including m = 2 and m = 3 hyperedges positive generalized hyperdegree correlations indicate for instance that nodes with many 2-body interactions have also many 3-body interactions and nodes with few 2-body interactions have also few 3-body interactions. On the contrary negative correlations of generalized hyperdegrees will imply that nodes with many 2-body interactions will participate in few 3-body interactions and viceversa. Secondly we might be interested in the overlap between hyperedges of different cardinality. This implies that in a hypergraph including m = 2 and m = 3 hyperedges, we might be interested to assess how many 2-body interactions connect nodes already connected in 3 body interactions.
In this work we will focus in particular in the effect of the correlations between generalized hyperdegrees on the robustness properties of hypergraphs. Indeed we notice that in the considered hypergraph ensemble hyperedges do not have a significant overlap. To show that that we define the total overlap O [m,n] between m-hyperedges and n-hyperedges with n > m as where Q [m] indicates the set m-tuples of nodes of the hypergraph and A [n] α = 1 if and only if α is a subset of nodes of an existing n-hyperedges, otherwise A [n] α = 0. The average overlap O [m,n] over the hypergraph ensemble with marginals given by Eq. (7) reads where p [m] α is given by Eq. (7) and where p [n,m] α is given by This implies that the average overlap O [m,n] is negligible for N 1 as it scales as (see analogous treatment for multilayer networks in [40]). Therefore the overlap of hyperedges is negligible in the sparse regime where the marginals are expressed by Eq. (7).

III. PERCOLATION ON RANDOM HYPERGRAPHS
Percolation in random hypergraphs can be treated directly by extending the ideas and concepts of percolation on factor graphs. Therefore, since the factor graph corresponding to a random hypegraph is locally tree-like we can write selfconsistent equations for the probabilityŜ that starting from a node and following a link we reach a factor node (hyperedge) in the giant component and for the probability S that starting from a factor node (hyperedge) and following a link of the factor graph we reach a node in the giant component. Assuming that each node is not initially damaged with probability p [N] and each hyperedge is not initially damaged with probability p [H] the self consistent equations for S andŜ read: A diagramatic representations of these two equations is shown in Figure 2.
The percolation problem is fully characterized by its order parameters given by the probability R of finding a node in the giant component and the probabilityR of finding a hyperedge in the giant component. In a random hypergraph, these order parameters can be expressed in terms ofŜ and S as These equations together with the self-consistent Eqs. (14) can be used to investigate the critical properties of percolation inferring the robustness of the random hypergraph. In particular we can impose p [H] = 1 (or p [N] = 1) and to characterize node percolation (or hyperedge percolation) where only hyperedges are randomly removed (or node percolation where only nodes are randomly removed). If the hypergraph only contains hyperedges of cardinality m = 2,( i.e. it reduces to a network) these two percolation problems reduce to link and node percolation respectively. In Figure 3 we show R versus p c . By imposing that the largest eigenvalue of the Jacobian matrix of Eqs. (14) is equal to one at S =Ŝ = 0, we find that p [H] c and p [N] c must satisfy for node percolation in which p [H] = 1 we obtain that the critical threshold p [N] c satisfies We note that Eqs. (17) and (18) fixing the hyperedge and the node percolation thresholds are invariant if we permute the distributions P(k) andP(m). This effect can be seen also in Figure 3 where it is evident that the hyperedge pecolation thresholds of two random hypegraphs with Poisson P(k) and PoissonP(m) is the same if for the first hypergraph k = 2 and m = 4 and for the second hypergraph k = 4 and m = 2.
Finally we note that the graph is formed only by mhyperedges, the distributionP(m ) reduces to a δ-function: and Eq. (16) reduces to: This last equations reduces to results obtained in Refs. [62,63].

A. General framework
We consider percolation on random multiplex hypergraphs, i.e. hypergraphs with given generalized hyperdegree sequences when nodes are not initially damaged with probability p [N] and hyperedges are not initially damaged with probability p [H] . In order to characterize percolation on these hypergraphs we consider their corresponding factor graphs. In particular we indicate withŜ m the probability that by following a link of a node in layer m we reach a m-factor node (mhyperedge) that belongs to the giant component. Moreover with S m we indicate the probability that following a link of a m-factor node (m-hyperedge) in layer m we reach a node in the giant component. Since the corresponding multiplex factor graph of the random multiplex hypergraph is locally tree-like the probabilitiesŜ m and S m can be find to satisfy the self-consistent equationŝ These self-consistent equations have a diagramatic interpretation as shown in Figure 4. In particularŜ m indicates the probability that a m-factor node (m-hyperedge) reached by following a link in layer m, is not iniitially damaged and it is connected at least to a node in the giant component. Instead S m indicates the probability that a node reached by following a link in layer m, is not initially damaged and it is connected at least to a factor node (hyperedge) -of any possible cardinalityin the giant component. The order parameters for percolation on a random multiplex hypergraph are given by the expected fraction of nodes R and the expected fraction of hyperedgesR in the giant component, given by: The Eqs. (22) together with the Eqs. (21) fully determine the percolation process in random multiplex hypergraphs and can be used to study the robustness of these structures as a function of the hyperdegree distribution P(k) and the distribution of cardinality of the hyperedgesP(m). In particular they can be used to investigate the effect that correlations between the hyperdegrees on different layers have on the robustness properties of the multiplex hypergraph. One fundamental measure for characterizing the robustness of a multiplex hypergraph with respect to another hypergraph is without any doubt the characterization of the percolation threshold. Indeed a smaller node (or hyperedge) percolation threshold implies that an hypergraph can display a giant component also when a larger fraction of nodes (or hyperedges) is removed.
The node and hyperedges critical thresholds p By simple analytical calculations it can be show that Λ is also the maximum eigenvalue of the matrix G of elements In the following sections we will predict the percolation threshold in important examples of random multiplex hypergraphs and we will characterize the role that correlations among hyperdegree of different layers have on the robustness properties of random multiplex hypergraphs.
B. Percolation threshold is some specific cases

Hypergraph with fixed cardinality of hyperedges
For a single layer multiplex hypergraph including only hyepredges of cardinality m, i.e. only including m-body interactions, we haveP In this case the matrix G reduces to a scalar G given by Therefore the percolation thresholds are obtained by imposing G = 1, giving It follows that in this simple case we recover the expression in Eq. (20) as we should.

Independent layers with Poisson generalized degree distribution
A more interesting case in which we can appreciate the multiplex structure of the problem is given by the case in which the hyperdegree distribution of each layer of the random multiplex hypergraphs is an independent Poisson distribution with layer-dependent average hyperdegree z m . In this case the joint hyperdegree distribution P(k) is given by with By using the well-known expression for the moments of Poisson distribution, we obtain that for this random multiplex hypergraph the matrix G has elements G mn given by Since the matrix elements G mn assume the same value of every element of a given row m the rank of G is equal to one, i.e. rank(G) = 1. This implies that the only non-zero eigenvalue Λ of G equals the trace of this matrix: By imposing that Λ = 1 we find that the critical thresholds  for each single layer of the multiplex hypergraph. Therefore the multiplex hypergraph is more robust than every of its layers taken in isolation.

Independent layers with power-law generalized degree distribution
Another interesting case of random multiplex hypergraph is the one formed by independent layers each one with powerlaw generalized hyperdegree distribution. In this case the joint hyperdegree distribution P(k) is given by with with γ m > 2 and c m indicating the normalization constant. For this random multiplex hypergraph the matrix G has elements Given that γ m > 2, we have that each layer is sparse, i.e. k n is finite at the limit N → ∞. However as soon as one layer is associated to a power-law exponent γ m ∈ (2, 3] the second moment k m (k m − 1) / k m diverges in the large network limit N → ∞. This implies that the trace of G diverges as well, indicating that the maximum eigenvalue diverges. It follows that as soon as one layer has a scale-free generalized hyperdegree distribution, i.e. as soon as for at least one layer m we have γ m ∈ (2, 3], then in the limit N → ∞. This implies that for standard percolation it is enough that one layer is scale-free to significantly improve the robustness of the random multiplex hypergraph.

C. Effect of correlations between generalized hyperdegrees
Random multiplex hypergraphs are characterized in general by non-trivial correlations between the hyperdegrees of the same nodes. In particular, given a random multiplex hypergraph we indicate with C mn the correlation between the hyperdegrees of the same node, connected to hyperedges of cardinality n and m respectively, i.e.
In a random multiplex hypergraph formed by two layers, this correlations can be modified by permutating the labels of the nodes in a given layer leaving the hyperdegree distribution of the two layers unchanged. In particular, it is possible to choose the permutation of the replica nodes in such a way that the correlations among the corresponding generalized degrees is maximized or minimized generating Maximally Positive Correlated Multiplex Hypergraphs and Maximally Negative Correlated Multiplex Hypergraph. This construction follows very closely the construction to build maximally positive and maximally negative correlated multiplex networks proposed in Ref. [15]. In particular the Maximally Positive Correlated Multiplex Hypergraph (MPCMH) can be obtained by ranking the generalized hyperdegrees of both layers in increasing order and identifying the label of the nodes with the same rank in both layer. On the contrary Maximally Negative Correlated Multiplex Hypegraph (MNCMH) can be obtained by ranking the generalized degree of one layer in increasing order and the one of the other layer in decreasing order, and by identifying the label of the nodes of the same rank. If the label of the nodes are assigned randomly we will obtain an Uncorrelated Multiplex Hypergraph (UMH). In order to assess the effect of correlations in the robustness of random multiplex hypergraph here we focus on a duplex hypergraph and we investigate the dependence of the percolations thresholds with the correlations coefficient between the generalized hyperdegrees of the two layers. In the considered case of a duplex hypergraph with two layers formed by hyperedges of cardinality m 1 and m 2 , the matrix G is given by wherem r = m r − 1 for r ∈ {1, 2} and where we have used the notation for r ∈ {1, 2}. The percolation threshold can be found by imposing that the maximum eigenvalue of the matrix G equals one, obtaining with We observe that for this percolation problem, the node percolation threshold p Since the product K 1 K 2 depends on the correlation coefficient C m 1 m 2 through Eq. (41) reveals that positive correlations increase the robustness of the multiplex hypergraph against random attack while negative correlations decrease the robustness of the multiplex hypergraph. In Figure 5 we show the effect of degree correlations on the robustness of random multiplex hypergraph by investigating separately node percolation and hyperedge percolation. We observe that the percolation threshold for node and hyperedge percolation are the same and are in perfect agreement with the analytical results indicating that maximally correlated multiplex hypergraphs have a lower percolation threshold than maximally negative correlated multiplex networks. ) calculated for the (MPCMH) and or for the (MNCMH). This implies that for large values of p the negative degree correlations enhance the robustness of the multiplex hypergraphs with respect to the positive correlations. In order to understand this phenomenon we note that close to the percolation threshold the robustness of the multiplex hypergraph is determined by the high degree nodes, that are less prone to damage in presence of positive correlations, leading to a smaller percolation threshold of MPCMH. On the contrary, for large values of p the robustness of the multiplex hypergraph, quantified by the fraction R of nodes in the giant component, is highly dependent on the low degree nodes. In particular the role of low degree nodes is more pronounced when in each layer there is a non neglibile number of isolated nodes. In presence of positive correlations among the generalized hyperdegrees, the number of nodes isolated in both layers or connected to a small number of hyperdeges (regardless of their size) is larger. As a consequence of this MNCMH have a larger fraction of nodes in the giant component than MPCMH, giving an intuition for explaining the fact that for large value of p the order parameter R becomes larger for MNCMH than for MPCMH in both for node and hyperdegree percolation in Fig. 5. This effect remains but it is strongly suppressed in absence of isolated nodes.  The topology of random multiplex hypergraph models allows us to explore a large variety of higher-order percolation problems. Higher-order percolation problems are characterized by illustrating cooperative phenomena where the probability that a node (or a factor node) is active depends on the presence of two or more active neighbours. These higherorder percolation problems have a highly non-trivial critical behavior and display hybrid discontinuous transitions, tricritical points, and they can even be characterized by more than one critical point as we will show in the following. Here we investigate and systematically characterize a large variety of higher-order percolation models that can be defined on multiplex hypergraphs. Inspired by the parallelism between multiplex hypergraphs and multiplex networks [6,37,39,64] we can define interlayer node interdependence in multiplex hypergraphs in which a node is active if it has at least a active neighbour on every layer of the multiplex hypergraph. This higher-order percolation is characterized by a hybrid discontinuous transition which can become a continuous transition at a tricritical point if partial interdependence is considered. However interlayer node interdependence is not the only interdependent model that can be defined on a multiplex hypergraph. In fact we can also consider interdependence associated to hyperedges, and assume that an hyperedge is active only when all its nodes are active. This highly non-trivial model display hybrid discontinuous transitions if the hyperedges are all involving more than two nodes. In presence of hyperedges of cardinality two (links) the transition can become continuous at a tricritical point in some cases. Note that this model is the percolation problem corresponding to the higher-order contagion problem proposed and studied in Refs. [11,26]. Another class of higher-order percolation problems is inspired by K-core percolation [35,36]. In the case of node K-core percolation a node is active if at least K of its hyperedges (of any given cardinality) are active, for hyperedge K-core percolation instead an hyperedges is active if at least K of its nodes are active. In either one of these last two models the transition is discontinuous as long as K > 2 and the distributions P(k) andP(m) have finite second moment. These different higher-order percolation problems are summarized in Figure 6.
B. Interdependent node percolation

General framework
In analogy to interlayer dependency on multilayer networks [6,37,39,64], we consider the interlayer dependency on hypergraphs. A node in the hypergraph is active when each of its replica nodes belongs at least to an active hyperedge, i.e. the node belongs to at least one active hyperedge for each possible value of the hyperedge cardinalities m. Therefore the probabilityŜ m that starting from a node we reach a m-factor node (m-hyperedege) that is active and the probability S m that starting from a m-factor node (m-hyperedge) we reach a node that is active follow the recursive equationŝ Moreover the order parametersR and R indicating the fraction of active hyperedges and the fraction of active nodes respec-tively are given bŷ As for interdependent percolation on pairwise multiplex networks, these equations lead to discontinuous (and hybrid) phase transitions. Let us indicate with Λ the maximum eigenvalue of the Jacobian matrix J of the system of Eqs. (44). The critical point of the discontinuous transition corresponding to non-zero order parameters R andR can be obtained by solving together with Eqs. (44) and Eqs. (45).

Independent layers
In order to reveal the mechanism responsible for the discontinuity of the transition, let us consider the model in the simple case in which the generalized degrees of a nodes are independent. In this case the joint distribution P(k) factorizes according to Eq. (27). In this limit the equations for S m and R can be simply written as where G 0,m (x) and G 1,m (x) indicates the generating functions By choosing the generalized degree distributions P m (k m ) to be Poisson and given by Eq. (28) these equations further simply as G 0,m (x) = G 1,m (x). Therefore we have R = S m = S for every possible value of R. Therefore in this simple limit the order parameter R = S obeys a single equations given by with the function h(S ) given by For any node-interdependent multiplex hypergraph with more than one layer this equation describes a discontinuous hybrid phase transition. Indeed the function h(S ) can display a minimum, when this minimum is achieved for S = S c such that h(S ) = 0 we observe the discontinuous phase transition (see Figure 7. Therefore the critical point can be found by solving We can consider the model in which p [N] = 1 or the model in which p [H] = 1. In both models the transition is hybrid.

Effect of correlations between generalized degrees
Interdependent multiplex hypergraphs display a higherorder percolation transition that is significantly affected by the correlations between generalized degrees of different layers. This phenomenon is the higher-order version of the corresponding phenomenon known to occur on pairwise multiplex networks [6,15]. By considering a duplex hypergraph with tunable correlations of the generalized degrees of the two layers we observe that MPCMH are more robust than MNCMH, i.e. positive correlations between generalized degrees of different layer increase the robustness of the multiplex hypergraph. This beneficial effect of positive correlations affects the critical threshold of the higher-order percolation model, which is lower for MPCMH than for MNCMH with the same hyperdegree distributions in each of the two layers (see Figure  8). Interestingly for interdepent multiplex networks the beneficial effect of positive correlations remains effective for every entity of the damage. In fact for this percolation problem, we have that also for large values of p [H] and p [N] the order parameter R for MPCMH remains always larger that the order parameter R for MNCMH. This phenomenology differs from the one observed for standard percolation (see Figure 5). The reason for this different behavior of interdependent percolation is simple: when most of the nodes and most of the links are not initially damaged, the fraction of active nodes is maximized for positive correlations. This remains true also in presence of isolated nodes. In fact negative correlations will imply the maximization of nodes which are isolated on at least one layer, and thus inactive, while positive correlations will minimize the number of nodes isolated in at least one layer. This simple explanation reveals why in Figure 8 the order param-eter R for MPCMP is always larger that the order parameter R for NPCMP, while we observe a crossing of the two curves for standard percolation (see Figure 5).

Partial interdependence
While node-interdependency always leads to discontinuous and hybrid transitions, if partial interdependence is taken into account it is possible to observe a change of behavior at a triciritical point separating a phase in which the percolation process displays discontinuous hybrid transitions from a phase in which the process displays continuous transitions. Partial interdependence has been introduced and investigated in detail for pairwise multiplex networks [6,39,64]. Here we extend this notion to multiplex hypegraphs highlighting the similarities and differences between the two models. By partial interdependence we mean that the interdependence is not always present between the replica nodes but replica nodes are interdependent only with probability r. Therefore for r = 1 we recover the node interdependent multiplex hypergraph studied in the previous paragraph and displaying a discontinuous hybrid transition, while for r = 0 we recover the standard percolation model studied in Sec. III displaying a continuous transition. Let us restrict our discussion here to the simple case of independent generalized degrees with join generalize degree distribution P(k) given by Eq. (27). In this case the equation forŜ m and the equation forR remains unchanged (given by the first of Eqs. (44) and (45), however the equations for S m and R change and are given by Interestingly due to the higher-order nature of the multiplex hypergraphs these equations cannot be reduced to a single equation in the case of Poisson layers with generalized degree distribution given by Eq. (28). However the phase diagram of the model can be investigated numerically. The phase diagram is characterized by a tricritical point separating a regime with r < r T for which we observe continuous transitions and a regime with r > r T in which we observe a discontinuous hybrid phase transition. Let us consider the case in which either nodes (p [H] = 1, p [N] = p) or hyperedges (p [N] = 1, p [H] = p) are randomly removed with probability 1 − p. In this case the tricritical point (r T , p T ) can be found numerically by solving the self-consistent equations forŜ m and S m together with where Λ is the largest eigenvalue of the Jacobian matrix of the equations determiningŜ m and S m (see Fig.9). Here p [N] is set equal to one, i.e. p [N] = 1.

General framework
In this section we introduce the higher-order interdependent hyperedge percolation model. In this model an hyperedge is active only if all its nodes are active as well and a node is active if at least one of its hyperedges is active. This model is here chosen because of its complementarity with the nodeinterdependence where a node is active if all its replica nodes are active, i.e. all its replica nodes belong to at least an active hyperedge. Interestingly the independent hyperedge percolation problem can be related to the model of higher-order social contagion proposed in Ref. [11] and investigated on random hypergraphs in Ref. [26]. Indeed in higher-order contagion model a node is infected if at least one of its m-hyperedges connects the new node to a set of m − 1 infected nodes. The inderdependent hyperedge percolation model and the higherorder contagion model can be mapped to each other. However there is a significant difference between the contagion model on hypergraphs and interdependent hyperedge percolation. While the contagion process admits a region of bistability in which the number of infected nodes can acquire either a larger or a smaller value depending on the initial conditions of the dynamics, in the corresponding region of the phase diagram, the interdependent hyperedge percolation does not display bistability. Indeed, although the self-consistent equation for the order parameter admits two solutions, the order parameter R always takes the value of the largest solution of the self-consistent equations. Broadly speaking percolation can be seen as an optimization problem in which one characterizes the maximum number of nodes that are connected under the condition imposed by the combinatorics of the process.
In the hyperedge interdepent percolation model the prob-abilityŜ m that starting from a random node we reach an mfactor node (m-hyperedge) which is active, and the probability S m that starting from a m-factor node (m-hyperedge) we reach a node that is active are given bŷ Moroever the order parameterR and R indicating the fraction of active hyperedges and active nodes respectively are given These equations differ with respect to the equation valid for standard percolation. In particular the equations forŜ m andR imply that an hyperedge can be active only if all its nodes are also active. Therefore we note that is the multiplex hypergraph contains only one layer and the layer capture only pairwise interactions, i.e. m 1 = m = 2 then this model reduces to standard percolation, however as long as the hypegraph contains hyperedges of cardinality m 2 the interdependent hyperedge percolation problem differs from standard percolation. In the following paragraphs we will investigate the nature of the percolation transition and the effect of correlations among generalized degrees observed for this model.

Independent layers
In is instructive to investigate the critical properties of hyperedge interdependence for a multiplex hypergraph with independent layers. In this case the Eqs. (54) and Eqs. (55) reduce tô where the generating functions G 0,m (x) and G 1,m, (x) are defined in Eq. (49). By considering Poisson layers with generalized degree distribution given by Eq. (28) we observe that S m = R = S for every value of m with S satisfying In the case of 2-layer multiplex hypergraph we obtain that S satisfies with Let us fix the expected number of hyperedges incident to a node, regardless of their cardinality, by imposing and let us investigate the nature of the interdependent hyperdge percolation transition as a function of z m 1 . Let us start with the specific example of having two layers with m 1 = 2 and m 2 = 3. If z 2 = z and z 3 = 0 the multiplex hypergraph reduces to a single network, and the transition is the standard percolation transition, which occurs at a critical point characterized by satisfying In the other extreme case in which z 2 = 0 and z 3 = z, the multiplex hypergraph reduces to a single layer hypergraphs including only 3-hyperedges. In this case the transtion is discontinuous and is obtained at a non zero value S = S c for which These are the two limiting cases of a region of the phase space in which we observe a continuous transition and of a region of phase space in which we observe a discontinuous transition. These two regions are separated by a tricritical point observed at the value of z m 1 = z T that satisfies For hyperedge interdependent percolation with p [N] = 1 we obtain the triciritcal point at For hyperedge interdependent percolation with p [H] = 1 we obtain the tricritical point at As we change the values of m 1 and m 2 characterizing the two layers of the duplex multiplex network different scenarios emerges. For m 1 > 2 and m 2 > 2, the transition is always discontinuous. Interestingly, as is shown in Figure 10(b), when m 1 = 2 and m 2 > 3 the hyperedge interdependent percolation can display not just one but also two percolation transition. The first transition describes the emergence of the generalized giant component and is continuous, the second transition indicates a discontinuity of the order parameter R from a nonzero value to another non-zero value. As far as we know this phenomenon has not been reported before, not even for the higher-order contagion model studied in Refs. [11,26] but can have an interesting interpretation in that context as a sudden activation of hyperedges of larger cardinality.

Effect of correlations
The general equations determining hyperedge interdependent percolation can be also used to study the effect of correlations between the generalized degrees of the replica nodes in different layers. In this case, regardless the nature of the phase transition, we observe that MPCHM display a transition threshold smaller than MNCMH, indicating that the system is able to sustain more damage. However for small entity of the damage, and in the extreme case in which the multiplex hypergraph is not damaged, the MPCMH have a smaller giant component than the MNCMH. This phenomenon is expected as it has the same explanation of the corresponding phenomena observed and discussed in Sec. III for the case of standard percolation (see Figure 11).

D. Node K-core percolation
In this section we propose the K-core node percolation on random multiplex hypergraphs. This model is a higher-order percolation process that generalizes K-core percolation of single pairwise networks to the multiplex hypergraphs. In K-core node percolation a node is active if has at least K active neighbours. In K-core node percolation defined on a multiplex hypergraph, a node is active if it belongs at least to K hyperedges regardless of their cardinality. The probabilityŜ m that starting from a node we reach an m-factor node (m-hyperedege) that is active and the probability S m that starting from a m-factor node (m-hyperedge) we reach a node that is active follow the recursive equationŝ where k indicates the sum over of k such m k m ≥ K.
Here B q (k) is given by where {q m } indicates the sum over of {q m } such that m q m = q.
The order parameters R andR expressing the fraction of nodes (R) and the fraction of hyperedges (R) in the node Kcore are given bŷ with D q given by It follows that the equations for node K-core percolation reduce to standard K-core percolation if the multiplex hypergraph if formed by a single layers encoding for hyperedges of cardinality m = 2 (i.e. links). For node K-core percolation, like for node K-core percolation on pairwise networks [35,36], we observe that the percolation transition is discontinuous and hybrid as long as K > 2 provided that the generalized degree distributions have finite second moment (see Figure 12).

E. Hyperedge K-core percolation
Hyperedge K-core percolation is here defined as a higherorder percolation process occurring on multiplex hypegraphs in which a hyperedge is active only if at least K (with K ≥ 2) nodes belonging to it are also active. In this case, the proba-bilityŜ m that starting from a node we reach a m-factor node (m-hyperedge) that is active and the probability S m that starting from a m-factor node (m-hyperedge) we reach a node that is active are given bŷ where B q (m) can be expressed aŝ Similarly we can define the order parameters R andR indicating the fraction of nodes and hyperedge that are active aŝ For hyperedge K-core percolation like for K-core percolation on pairwise networks [35,36], we observe that the percolation transition is discontinuous and hybrid as long as K > 2 provided that the distributionP(m) has finite second moment (see Figure 13).

VI. CONCLUSIONS
In this paper we have provided a comprehensive framework to study standard and higher-order percolation on random multiplex hypergraphs. Random multiplex hypergraphs are a natural generalization of random hypergraphs where the hyperedges of different cardinality are associated to different layers of the multiplex. This modelling framework is very comprehensive and is here used to investigate the rich interplay between the topology of hypergraphs and the properties of standard and higher-order percolation defined on these structures. We reveal how interlayer correlations among the generalized degree of replica nodes can affect the critical properties of standard percolation. In particular we show that close to the percolation transition positive correlations enhance the robustness of multiplex hypergraphs while when the initial damage is minor negative correlations can be beneficial to network robustness. We show how the multilayer nature of multiplex hypergraphs can be exploited to define a number of higher-order percolation processes. In particular we propose two models generalizing interdependent percolation in multiplex networks and contagion model in hypergraphs (the interdependent node and the interdependent hyperedge percolation) and two models generalizing K-core percolation to hypergraphs (the node K-core and hyperedge K-core percolation). These models are here shown to display a rich phenomenology including discontinuous hybrid phase transitions, tricritical points, and multiplex phase transitions together with non-trivial effects due to the interlayer correlations among the generalized degrees.
Although our aim is to provide a comprehensive view of the possible higher-order percolation processes on random multiplex hypergraphs we are aware that the processes investigated in this work are not exhausitive of the many relevant percolation processes that can be defined on these structures. We hope that this work can generate further interest in the interplay between the structure of higher-order networks and their dynamics and that the revealed properties of percolation on multiplex hypegraphs can open new insights also for the study of other dynamical processes such as epidemic spreading and social contagion.