Interaction Measures, Partition Lattices and Kernel Tests for High-Order Interactions
Authors: Zhaolu Liu
Robert Peach
Pedro A. M. Mediano
Mauricio Barahona
Venue: NeurIPS 2023 poster (2023)
Keywords: High-order interactions; Lattice theory; Kernel tests

Abstract: Models that rely solely on pairwise relationships often fail to capture the complete statistical structure of the complex multivariate data found in diverse domains, such as socio-economic, ecological, or biomedical systems. Non-trivial dependencies between groups of more than two variables can play a significant role in the analysis and modelling of such systems, yet extracting such high-order interactions from data remains challenging. Here, we introduce a hierarchy of $d$-order ($d \geq 2$) interaction measures, increasingly inclusive of possible factorisations of the joint probability distribution, and define non-parametric, kernel-based tests to establish systematically the statistical significance of $d$-order interactions. We also establish mathematical links with lattice theory, which elucidate the derivation of the interaction measures and their composite permutation tests; clarify the connection of simplicial complexes with kernel matrix centring; and provide a means to enhance computational efficiency. We illustrate our results numerically with validations on synthetic data, and through an application to neuroimaging data.

Interaction Measures, Partition Lattices and Kernel Tests for High-Order Interactions
Zhaolu Liu'
Robert L, Peach
Pedro A.M. Mediano
Mauricio Barahona"
Department of Mathematics, Imperial College London, United Kingdom 2 DDpartment of Neurology, University Hospital Würzburg, Germany Department Brain Sciences, Imperial College London, United Kingdom ' Department of Computing. Imperial College London, United Kingdom
Abstract
Models that rely solely on pairwise relationships often fail to capture the complete statistical structure of the complex multivariate data found in diverse domains, such as socio-economic, ecological, or biomedical systems. Non-trivial dependen- cies between groups of more than two variables can play significant role in the analysis and modelling of such systems, yet extracting such high-order interac- tions from data remains challenging. Here, introduce hierarchy of d-order interaction measures, increasingly inclusive of possible factorisations of the joint probability distribution, and define non-parametric, kernel-based tests establish systematically the statistical significance of d-order interactions. We also establish mathematical links with lattice theory, which elucidate the derivation of the inter- action measures and their composite permutation tests; clarify the connection of simplicial complexes with kernel matrix centring: and provide means enhance computational efficiency. We illustrate our results numerically with validations on synthetic data, and through an application to neuroimaging data.
Introduction
There is increasing evidence that pairwise relationships are insufficient to model many real world systems 2
3 The relevance of high-order interactions has been emphasised in many contexts, as relationships within social [], ecological [S] 6), and biological systems [7] frequently involve groups of three or more agents, beyond pairwise associations. Such high-order interactions can be neither trivially represented by linear combination of dyadic relationships, as the presence of high-order interactions can significantly impact the dynamics on networked systems
10
13 14 15).
nor can they be simply detected by joint independence tests which inherently ignore other possible factorisations of the joint distribution.
Direct measurements of group interactions are seldom available, leading to their omission, partial representation as projected pairwise interactions [16] 17 [8] limited recovery through inference
methods [19 20] Often, only indirect measurements of underlying interactions in real-world complex systems are available, the form of iid or time-series data. Methods that learn high-order interactions from time-series data have been developed [I]] 21), however, their reliance on heuristic measures makes their interpretation difficult. Alterative methods have been developed with strong foundations
in statistics [22 and information theory [3 24 25), however, as we show later, these are often based
Corresponding author: m.barahona imperiallac..u
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

d 2
d 3
d 4
additional Steeiberrg tarms
dierrbuttion
debiburiis
partill
compeete Gotion
complett
Figure I: Factorisations of joint distributions Fi...d for d 2,3, The black dots indicate the marginal distributions of the singletons. The line, triangular and square shapes represent the joint distribution of two, three and four variables respectively. Different factorisations are presented as partitions of the d variables ordered by increasing cardinality from top to bottom, so that all the factorisations with the same number of independent blocks appear the same level. Joint independence considers only the top and bottom levels for each d. whilst the Lancaster interaction considers all terms except those in the shaded region. The Streitberg interaction considers all partitions. A whereas for d = 3, we have AjP + A A
Hence, for d = 2, we have 7P A P-
and for d 2 4, we have P + A P + Ag PP.
on an incomplete set of interactions and thus fail to capture all possible factorisations of the joint probability distribution [26 27].
Kernel-based hypothesis tests provide a non-parametric, statistically robust framework for detecting relationships between variables from observational data. Such tests have been implemented for pairwise 28 29 30] [1] and d-variate joint independence [32] 33 J. and proven to be effective for
non-trivial dependencies such as in the 3-way Lancaster test 34 35 However, the best of our knowledge, cases where the number of variables exceeds 3 are still unexplored.
Here, extend the capability for detection of high- nigh-order interactions by introducing family of tests based on factorisations of the joint probability distribution that generalise systematically to any order d. At the head of this family, introduce the Streitberg interaction test, which captures all factorisations of the joint distribution
order d. We further show that, despite the fact that the naive extension d 2 4 of the Lancaster interaction excludes some of the factorisations of the joint distribution, not all is lost, and we detail the conclusions that can be drawn from rejecting the d-order Lancaster interaction. Furthermore, foundation for detecting high-order interactions; (ii) interaction measures can be systemically derived from partition lattices; (iūi) the corresponding Hilbert-Schmidt norm for kernel embeddings can naturally be extracted from the product lattice; (iv) the composite permutation tests can be performed with regard the second level of the lattice; and (v) the lattice formulation allows us to propose a generalised interaction measure that can be used to test whether given factorisation can be factorised further. Despite the inescapable combinatorial nature of testing high-order interactions, we offer approaches to reduce the computational complexity of our d-order interaction tests informed by our links with lattice theory. Finally, we present numerical validations of our tests on synthetic data before applying them to real-world neuroimaging data.
show that: (i) lattice theory provides coherent theoretical
2 Interaction Measures
The most basic form of interaction between variables occurs between two variables, X and X and is often characterized by the lack of pairwise independence. i.e., the difference between the joint distribution xxx*,x and the product of the marginal distributions Fx1 and Fxxa (see Fig D 2).
Let us consider d random variables {X' X²_ ......... Xd) with joint distribution Fxx..... (for readability, hereafter we use this subscript notation). The d variables are jointly independent if and only if Pi can be factorised as the product of the univariate marginal distributions. Hence we can introduce an interaction measure that vanishes only when the variables are jointly independent:
i
d
A/P = Pl...
II P.
When d > 2, however, the criterion for joint independence fails to capture high-order interactions, since additional factorisations must be considered. For instance, if 122 could be factorised as

---

P1P23 (see Fig. 1 d d= 3). joint independence would be rejected, yet there is no 3-way interaction. Instead, we seek way of identifying high-order interactions for variables where all lower order independencies can be rejected.
One approach to address the shortcomings of joint independence in identifying high-order interactions for d > 2 variables is to use a signed measure called the Lancaster interaction [34]. The Lancaster
interaction for d = 3 is defined as A 123 1 F22l113 32'1''1 12 + 2P 2P A J. vanishes if IP 22 can be factorised in any way, and has been implemented as a kernel-based test statistic to
identify factorisrie joint distributions 33 35 Lancaster also generalised this measure the
multivariate case with d variables:
A P
(IP) - P:),
==1
where PPP ....P.. PP Pij.-.k The Lancaster interaction vanishes when the joint distribution can be factorised into jointly independent subvectors for d = 3. but it fails when d 2 [26]. Specifically, A P does not vanish if P1244 factorises into PizIPs4, F1 12 Despite this, the Lancaster interaction is not entirely uninformative for d 2 4, and in Section 4.1 we examine the necessary
conditions for it to vanish.
The desired vanishing property can be achieved by using the Streitberg interaction 26
more general
interaction measure defined using partitions. Let 7) be the set of random variables (X X² Xd), and let (I) denote the set of all partitions of D, where a partition π is a collection of nonempty, pairwise disjoint subsets (blocks) bj S D that cover D. Then, the Streitberg interaction measure is defined as
P
Σ (()) ()1
TELI(D)
Here, |/|| denotes the cardinality of the partition T, and P 1 Po, is the corresponding factorisation with respect to b1 [1]62]... 1b,, It has been proven that A 0 if the joint distribution can be factorised in any way, although the converse is not true in general 33
3 Partition Lattices
The expressions of the interaction measures outlined in the previous section be systematically generated from partition lattices. Interestingly, this formulation also allows us to establish further theoretical links with simplicial complexes.
Notation. A partially ordered set defined on a set S with a binary relation < is a lattice C if for any
a L there exists a greatest lower bound (meet) a A and least upper bound (join) V π (36] and Ó. Let A--) denote real-valued
We denote the maximum and minimum element of a lattice function defined on S, then the sum function f(-) is analogous to the integration of A(-) over the interval (0, r.. where æ E L. Then A(- can be expressed as the inverse operation of f()):
f(T) Σ (((am))
A(T) = Σ M(o, x) f(a),
(4)
where the partial order is encoded by the Zeta matrix, with <(a,T) = of < and otherwise. Its inverse is the Möbius matrix with elements pla, σ) (for details see Section Dj. which can be obtained explicitly 37.
Lattices and Interaction Measures. The construction of the interaction measures in Section 2 is closely related to the partition lattice [26] 37). The partition lattice is defined on the set IT(D) with ordering given by the notion of partition refinement. A partition denoted as < æ, if every block of a is fully contained within a block of π. The lattice structure thus allows us define the least upper bound a and the greatest lower bound π between any two partitions o and π. Clearly, the maximum element of the partition lattice, î, corresponds to the joint distribution, whereas the minimum element, Ó corresponds to the complete factorisation (Fig. 1
said refine another partition x,
Within this formalism, the interaction measures are obtained when the sum function f() is the probability distribution function. It then follows that the interaction measures are obtained
3

---

from the inverse operation. Indeed, the Streitberg interaction (3) is given by the Möbius inversion defined over the complete partition lattice. In contrast, the Lancaster interaction is obtained when the inversion is defined over the subset of partitions that have at most one non-singleton block, which we denote as the Lancaster sublattice. In other words, the sum function associated with the Lancaster interaction considers fewer interactions than the Streitberg interaction and its Möbius inversion has correspondingly fewer terms, as seen in 2p. Note that for d = 3, the full (Streitberg) lattice and the Lancaster sublattice are the same, hence the interaction measures coincide. For d 2 4. however, the Streitberg and Lancaster lattices, and consequently the corresponding interaction measures, are different, shown by the extra partitions with two non-singleton blocks (shaded region) in Figure [ Note also that joint independence considers trivial sublattice with only two elements: Ô and (for any d). The Möbius inversion on this sublattice leads to (T). which only vanishes for complete factorisations.
Links Simplicial Complexes. A popular representation of high-order systems in the literature is through simplicial complexes (38]. Importantly, the simplicial complex construction can also be understood in terms of partition lattices. In particular, the elements in (d simplex have inclusion ordering and thus form a subset lattice, e.g.. {X) {
} is subset of { xx2, X The
subset lattice has been utilised in the understanding of information geometry [3] and solving tensor balancing [40] It can be shown that the deatomised sublattice (with removal of singletons) is
isomorphic to the Laneaster lattice 41). Furthermore, the boundary matrices of a (d simsmp appear as block matrices both in the Zeta matrix and the (inverse) Möbius matrix of the d-order partition lattice. These are explored further in Section in the SI.
Note that the non-singleton partitions, i.e., those that do not belong to the Lancaster sublattice, are in the set of simplicial complexes. Hence these blocks and their respective refinements cannol be expressed in terms of boundary matrices. Therefore the full partition lattice and its resulting Streitberg interaction provides more information compared with measures originating from sublattices, and in particular the Lancaster lattice and its related simplicial complex construction [23] 24]. Whilst simplicial complexes, and similarly the Lancaster interaction, can recursively constructed from lower order elements, this is not possible for the Streitherg interaction due to the partitions without singletons, which has implications on computational efficiency (Section 5).
4 Kernel Interaction Tests
The interaction measures in Section can be utilised as statistics in non-parametric tests when embedded into reproducing kernel Hilbert spaces (RKHS). Given symmetric, positive definite function le : d' x X + R. there is an associative RKHS H' with the reproducing kernel property. For X E X
denote (0(() the canonical map of ks (X .). The kernel mean embedding of Satisfies is characteristic, the mean embedding is injective and the norm of the signed measure is zero if and only if the measure is zero itself [28] 42 These properties enable us to create meaningful non-parametric tests by computing the kernel mean embedding of desired interaction measures.
Pxs, Me Mzx. satisfies Exx)(X') (f, uz xi /1s where HS stands for Hilbbett-Schmitt If the kernel
In this section, we extend interaction measures for random variables the d-variate case, and formulate family of interaction tests including the Streitberg interaction, the Lancaster interaction, joint independence and a generalised interaction.
For proofs of the propositions, please see Appendix
4.1 Lancaster Interaction
Let us first consider the Lancaster interaction. Although it does not necessarily vanish when d' 2 4 for all types of factorisations due the lack of certain partitions, here we find the necessary conditions for the Lancaster interaction to vanish.
Proposition 1. If the joint distribution Pi.--d can be factorised into Px. where ny are partitions with at least one singleton, then A O.
Remark: Note that IPr,, is broader set of partitions compared to the set of partitions in the Lancaster lattice e-g., for d = 5. IP12IPasPs P3 satisfies Proposition| [] yet it is not constituent of the Lancaster lattice, which only consists of partitions with at most one non-singleton.

---

We now define the Mixed Central Moment Operator as the kernel embedding of the Lancaster interaction, which follows immediately from the expansion of 2 (27 Definition (Mixed Central Moment Operator).
d
Ma=
E)
+ Tep0
(-1)) E
II
00 II
[0)
(5)
where he are the set of partitions with most one non-singleton (i.e., those that belong to the Lancaster lattice), are the singletons, and j runs over variables in the non-singleton block.
Proposition 2. By rearranging, the Mixed Central Moment Operator can be simplified to:
Ma =E=
II o' -
(6)
Remark: This simplification transparently re-expresses the operator as central moment, instead of the complex sum of moments (5). The simplification eliminates all partitions except 1.
We now define the entries in matrix Kib kt (x szf)) for vid samples xa and æf where S a, b < n
and K = НK*H H where H = I - 111 is the centring matrix. The norm of the embedding above can serve test statistic. The estimator of the test statistic derived as the V-statistic: Proposition 3 (Lancaster interaction estimator.).
d
Malless
2b>
(7)
n
4.2 Streitherg Interaction
Similarly we define the Mixed Cumulant operator in terms of kernel embeddings from the Streitberg interaction as:
Definition 2 (Mixed Cumulant operator).
K.a
Σ (-1- TEI((D)
(8)
bE
iEb
where is element in block of partition in partition lattice II(D).
This operator be simplified in similar way:
Proposition 4.
K Σ (1T (|(-)-1)((1))--))
E
II LiEb
06 E[d']
(9)
E.E(D
DEF,
where n, are the partitions with no singletons in partition lattice II(D).
Remark:
Note that partitions with singletons are eliminated after centring, and the remaining partitions without singletons form an upper semi-lattice. The relationship between the Mixed Cumulant operator and the Mixed Central Moment operator is given by:
The Mixed Cumulant operator is equal to the sum of Mixed Central Moment Operator products associated with partitions with no singletons denoted as n. (27]:
Lemma 1.
K-d
Σ (1T.l (1m_1-1)(--1)) 1)1(-1)
M
(10)
N. E.I(I
DER,
Remark:
When d 2,3, Ka and Ma are identical, which further confirms the moment-cummlant relationship.
Remark: From the partition lattice (0. C ग. C म, where म. the set of partitions the Lancaster lattice and æ is the full set the Streitberg interaction. When d = 2, we have {Ó. i} = Me = म, and
when d = 3, we have (Ō, 1} C क n.
The estimation using V-statistic for the norm of the Streitberg interaction is given by:
5

---

Proposition 5 (Streitberg interaction estimator),
Σ
(IT=| - 1))(1==((== - 1)((-1)(sl+t))
1 n|al++/51
ΣΣЛ (be} 16
K b
(11)
|=
where 1 S bi, and [(b,]] stl [8 |T'')- i.e. indices b, by if i.j are in the same
block in na similarly for b,. The sets of indices {bi}, {6,} are disjoint. An example for d = 4 is given in the ST.
Norms and the Product Lattice. The test statistics of the interaction measures, i.e., the norm the interaction operators, can also be directly obtained from the Möbius inversion of the Cartesian product of the upper semilattice, thus avoiding the need to compute the square of the operator. The Cartesian product of two lattices with product ordering is also a lattice 43]. The elements in the Cartesian product are analogous to the inner products of kernel embeddings in the estimators. Therefore, the coefficients can be computed directly from the product lattice, as shown in Appendix
4.3 Description of the Statistical Tests
The Streitberg interaction test involves rejecting the null hypothesis that the joint distribution can be factorised in some way, which occurs when A o O. Similarly, for the Lancaster interaction test, we reject the null hypothesis that the joint distribution can be factorised into partitions with at least one singleton by checking A P * O. Both null hypotheses consist of multiple sub- sub-hypotheses, and
an example for d = 3 is discussed in Ref. 35). Importantly, the number of sub-hypotheses to be tested is related factorisations with only two blocks, i.e., those that form the second level from the top in the partition lattice (see Figure [ In the case of rejecting some but not all the hobthhyttsesss, please see Appendix for the discussion.
To test each sub-hypothesis II.... PoPv. for efficiency we fix the observations of the variables in the largest block, and permute the observations of the remaining variables in the other block to induce the independence structure in the null distribution. We then use the Monte-Caria approximation to compute the p-value [32], and apply simple correction 335 generalised to d variables to correct for multiple hypothesis testing. Finally, we reject the composite null hypothesis if all corresponding sub-hypotheses are rejected. To achieve this, we employ the permutation test outlined in Algorithm
Algorithm 1 Permutation test for the interaction measures
test-statistic || [| (K.....K
for i E {m]} (here л = bjD') do
initialise empty P-dimenniinal vector T
for 1: do
for j = 1: d do
if argınin(<<]) then zE (b.b')
Ki K2 (c) (kernel matrices after random permutations the observations of X
Tpp |1
Rd)
tmp + {p E {1. PH Tpp 2 test-statistic)
pval <- (tmp + 1) 1) approximation 22
if pval > then
reject < 0 ( ssimplee correction' ' 35
break
reject +
return reject
By rejecting the null hypothesis of the Streitberg interaction test, it follows that we reject the null hypothesis of the Lancaster interaction test, as well as the null hypothesis of joint independence. The set of sub-hypotheses in joint independence is a subset of the set in the Lancaster interaction test, which is in turn subset of the set in the Streitherg interaction test. The number of sub-hypotheses for an interaction measure is equal the number of elements the second level of the respective lattice.

---

For Streitberg interaction, this number is (2 1) as follows from the full partition lattice; for Lancaster interaction, this number is d due to the reduced Lancaster lattice: for joint independence the number of sub-hypotheses stays fixed as the corresponding 2-element lattice always contains precisely one element, i.e., corresponding to P1 P. Pa.
Remark: Although our focus is composite interaction tests unveil high-order interactions by leveraging the vanishing of Streitberg and/or Lancaster interactions, these measures are also valuable test other lower-order independence hypotheses wherein factorisations can also lead to vanishing interaction measures. For example, both Lancaster and Streitberg statistics can be implemented to test for joint independence (as shown in Section 6). and for marginal independence.
4.4 Generalised Interaction Measure
The measures described so far are only able to handle interactions within one group of variables. But what if are interested in knowing whether a factorisation (e.g., P12P31) factorises further? To answer this, we can simply formulate an interaction measure based on a partition [44]. Because an interval of lattice L is subset of the form (o, x = {z interval lattice can be used to produce generalised interactions for multiple groups of variables:
л and is also lattice [4]], the
Definition 3 (General interaction operator).
A 5 P = Σ
m(a, T.P., = II A
(12)
TST
DET,
where na are the partitions with no singletons, m are the Möbius coefficients, and b are the blocks in the partition T.
Here, we only consider blocks of size at least 2, since a factorisation cannot not be defined for singletons. It then clear that A IP 0 if any A&P is zero, and the rejection statement on each group will be dependent on their own vanishing conditions. This provides a general framework to detect interactions within blocks of variables. The corresponding kernel-based tests can be constructed by estimating the norm of A5 P. In fact, the generalised Streitberg interaction measure can be expressed in terms of ordinary Streitberg interaction measures defined in Section 4.2 44
A S
P
Σ I eVR,=ibEc DEo
(13)
5 Computational Considerations
The time complexity of computing the Lancaster interaction estimator is O(d))), where d is the order of interaction and is the number of samples. This follows immediately after centring the kernel matrices, since is the only partition with no singletons. In contrast, the Streitberg interaction has a larger number of partitions, given by the Bell number, B& [45). However, by simplifying the computation through centring, can reduce this the number of partitions with no singletons, given by Fa [46]. While this number is still combinatorial, it is significantly smaller compared to Ba (see Table in Appendix G).
The resulting estimator of the Streitberg interaction criterion then consists of Fa inner products. The summation indices of the inner products come from both partitions and are most d (when d is even) or d- (when d is odd). Whilst the summations can be contracted, the choice of contraction ordering plays crucial role in the time complexity [47 48 49 By computing the summations related to one partition in parallel, since the index sets from n. and T are disjoint, the overall complexity can be reduced to On "0m= (see Table 2 in Appendix 5). For example, the inner product between
the kernel embeddings for P1234 and P12I'' 34 is n Σ. E, Σx K! K? K3 K& but an improved contraction ordering is + E, KIJKII,,
J
The time complexity of the Streitberg interaction estimator can be further reduced in a recursive setting. Specifically, if n. V n', < i, then the associative inner product can expressed as a product of independent sums that can be found in lower order interactions. For example, (14 (araFaa) can be decomposed as (MP a (2rap) (UTay (Paa)) which would have already been computed for the pairwise interactions of {X = and (X3, x.) This means that if one has already computed
7

---

a
C
Factoriaatin: P13P3445
Lancact:
aMsic
atAc
Lancsee
1.3
mekito uesiChvedes diach
Statler Sossllerghutt
0.2
25
03
0.2
1.2
Fleection
Inteactie arsigh
Interacc puppertinn
Uiterct properion
Figure 2: High-order tests on synthetic data. (a) dHSIC unable to detect the marginal P1 P2243 factorisation. (b) Only the Streitberg test is able to capture the factorisation without singletons 5-way XOR example, Lancaster and Streitherg tests achieve higher power than dHSIC for the same number of samples (n=80). Only by increasing the number samples to n=1000 does dHSIC display comparable power. (d) The modified dHSIC test requires substantially more samples than the Streitberg test display comparable power when testing all sub-hypotheses the interaction in the XOR example.
P12 (c) Testing joint independence
the inner products of d - 2 variables (incremental by 2 since the smallest block is at least size 2), then the only remaining inner products to compute are between those embeddings associated with partitions я, and nT' such that M, V n. i. other words, the computation of inner products for high-order interactions can be reduced to the computation of lower order interactions, long the corresponding partitions can be joined to form a partition of all d variables. For further discussions on computational complexity, see Appendix G
6 Experiments
We first validate our family of d-order interaction tests (in particular, d 5) synthetic datasets with ground truth interactions, and then investigate their application to a neuroimaging dataset. Unless stated otherwise, the significance level is set to = 0.05, sample size to n 80, and we use Gaussian kernels with the median heuristic as the bandwidth.
Synthetic Experiments: Multivariate Gaussians. We first compare results from apply- ing the dHSIC test (developed solely to detect joint independence) [32. and the Lancaster and Streitberg interaction
five-variable multivariate Gaussian distributions. N(n, = [0, 0, 0,0,0.00 and Σn (3,11 f,, interaction strength between variables and consequently P12315 = 12115 For the second example, have the same μ and ΣI = [[L, 3, 0), [8, 1, 1,0, O, 0], [0, O, 1, 8, 8), [0, 8, 1, 3), [0, 0, 8, 1]]
with mean fa and covariance Σ. For the first example,
[[1,0, O, 0, 0): [0, ,,, $) [0, 3, 3, 33 10,3,,3 9, 3), [0, B,
where 0 9 < 1 defines the
and hence P12345 = 12221
Figure a-b shows that dHSIC fails capture the partial factorisation of the joint distribution in both cases. Whilst the Lancaster interaction test is able to detect the factorisation in the first experiment (Figure 2(a)), it fails on the second experiment where the factorisation contains singletons (FigureD) b)). The Streitberg interaction test the only one that detects both factorisations.
Synthetic Experiments: 5-Way XOR Gate. Next, we investigate the power of the tests using
data set constructed using XOR gate. We generate samples of V, X, Y, i...d. U(0, 4). Z. = (V., + W + X,+ + Y,)) mod 4 and Zi+1m ~ (0, 4)) (where the samples [ 1: ni] act as noise). We then gradually increase the interaction proportion, 0 < i/n < 1. By construction, this dataset does not contain pairwise, 3-way or 4-way interactions, and the 5-way interaction becomes increasingly easier detect i/n approaches 1.
Given that the joint distribution this example cannot be factorised, all tests should reject joint independence. We show in Figure that the rejection rates of joint independence for Streitberg and Lancaster approach one as the interaction level increases for a small number samples (n 80), whereas dHSIC displays low power for the same number of samples. Only by increasing the number
8

---

diander
2-way
HD) petter Lage
0.-
0.0
Figure 3: Concentration of high-order interactions in brain resting state networks. The per- centage of rejected null hypotheses for d-order (d = 2, 3, 4, 5) interaction tests, when the d regions are selected within each of seven brain resting state networks (RSNs) or drawn at random across the brain. For each order d, we carry out Fisher-exact tests for repreeepresnttation for each RSN vs. random. All Fisher tests are sigmificant (p-value<< 0.05) except where shown as (n.s.). Abbreviations for RSNs: DMN: default mode network, SOM: somatomotor, VIS: visual, SAL: salience, DAN: dorsal attention network, FPN: frontoparietal network, LIM: limbic.
of samples to 1000 do we observe comparable performance for dHSIC. To test for 5-way interactions the XOR data set, we omit the Lancaster test which cannot detect factorisations into non-singletons. We thus compare the Streitberg interaction test with modified dHSIC [34 35 that tests each individual sub-hypothesis (e.g.. when testing (X,Y) 1 (Z, W) sufficient to treat Y) and (Z, W) single variables and perform joint independence tests for two variables). We find again that the modified dHSIC requires a much larger sample size to achieve comparable performance that of the Streitberg interaction test (Fig. dd)
To investigate how the test power degrades when the sample size is decreased, we performed further experiments in Fig 8 in Appendix F In all cases, we find that the performance of dHSIC degrades faster than that of Lancaster and Streitberg, i.c., the dHSIC null rejection rate decays substantially faster than that of Lancaster and Streitberg the sample size decreases.
Neuroimaging Dataset. As a proof of concept of applications to real data, we apply the Streitberg interaction test to detect high-order interactions in brain activity data. The dataset consists of resting- state fMRI data from 50 unrelated subjects part of the Human Connectome Project [50] [1]. For each subject, the data consists of 100 time series capturing the activity of each of the regions in the Schaefer brain atlas [52]. Importantly, these 100 regions be divided in groups called "Resting State Networks' (RSN) [53] Experimental evidence has shown that regions within the same hypothesised that high-order interactions would be more
RSN perform similar functions 52
probable within RSNs.
To avoid dealing with issues of non -on-staiionarity in the time-series data, we first take temporal averages obtain the mean activity per region and subject (see Appendix for preprocessing steps). We then perform the Streitberg interaction test to sets of brain regions that are either taken from within the same RSN or sampled at random from the whole brain. In each case, we sample 500 sets of regions or take all possible combinations, whichever is lowest.
The results in Figure 3 align with our general hypothesis: 2-, 3-, 4- and 5-way interactions are significantly more common among regions belonging to the same RSN as compared to random sets of regions, confirming that our interaction test can successfully detect high-order interactions in real-world data. Notably, our results show a larger percentage of 2-way interactions in the SOM and VIS compared regions such as the FPN or DAN, suggesting that there is increased redundancy in the system relative to these other regions. The SOM and VIS are structurally coupled, modular sensorimotor processing regions, that benefit from information redundancy to increase robustness. Moreover, the proportion of high-order interactions (3-, 4- and 5-way) to 2-way interactions is larger in SOM and VIS compared to FPN or DAN, potentially reflecting the differing balance of redundancy and synergy in these brain regions (54). Additional analyses are shown in Appendix H. but the

---

neuroscience implications of these results, including observed variation across regions, will be the object of future work.
7 Discussion
In this work, we have established formal connections between lattice theory and a family of high-order interaction measures based on factorisations of the joint probability distribution. The link lattice theory facilitates the derivation and interpretation of the measures; highlights the connection to simplicial complexes and other sublattices; and helps formulate the associated permutation tests more efficiently. We have shown empirically that the Lancaster test statistic is good substitute for dHSIC test for joint and marginal independence, whilst the Streitberg test statistic is able to better capture all factorisations of the joint distribution. Our work offers a rigorous and systematic approach for the reconstruction of high-order systems, but has some limitations. The computation of the Streitberg interaction estimator can be expensive and although we have proposed strategies to bring down the time complexity, the number of terms remains combinatorial (see Appendix G discussions on practical limitations). In addition, our theoretical results rely on iid data, a strong assumption in many real-world applications. This leaves several open directions. investigate whether the null distribution of Ag P can be jointly approximated without the need for sub-hypote hoto-hyotthees, thus reducing computational cost, how it can be accurately approximated when the data has temporal dependence; and what role high-order interactions may play causal discovery.
particular, future work will
Acknowledgements
MB acknowledges support by the EPSRC under grant EP/NO14529/1 funding the EPSRC Centre for Mathematics of Precision Healtheare at Imperial, and by the Nuffield Foundation under the project "The Future of Work and Well-being: The Pissarides Review". RP acknowledges funding from the Deutsche Forschungsgemeinscaatf (DFG, German Research Foundation) Project-II> 424778381-TRB 295.
For the experiments in Section 6 data were provided by the Human Connectome Project, WU-Minn Consortium (Principal Investigators: David Van Essen and Kamil Ugurbil; 1U54MH091657) funded by the 16 NIH Institutes and Centres that support the NIH Blueprint for Neuroscience Research; and by the McDonnell centre for Systems Neuroscience at Washington University.
The authors would like to thank Jianxiong Sun, Xing Liu and Paul Expert for valuable discussions, Asem Alaa for help in maintaining the GitHub repository, and Andrea I. Luppi for assistance with the preprocessing of neuroimaging data.
References
[1] Federico Battiston, Enrico Amico, Alain Barrat, Ginestra Bianconi, Guilherme Ferraz de Arruda, Benedetta Franceschiello, Iacopo lacopini, Sonia Kefi, Vito Latora, Yamir Moreno, et al. The physics of higher-order interactions in complex systems. Nature Physics, 17100):0033-109, 2021.
[2] Federico Battiston, Giulia Cencetti, lacopo lacopini, Vito Latora, Maxime Lucas, Alice Patania, Jean-Gabriel Young, and Giovanni Petri. Networks beyond pairwise interactions: structure and dynamics. Physics Reports, 874:1-92, 2020.
[3] S. Boccaletti, P. De Lellis, C.I. del Genio, K. Alfaro-Bitter, R. Criado, S. Jalan, and M. Ro- mance. The structure and dynamics of networks with higher order interactions. Physics Reports, 1018:1-64, 2023.
[4] Alice Patania, Giovanni Petri, and Francesco Vaccarino. The shape of collaborations. EPJ Data Science, 6:1-16, 2017.
[5] Peter A Abrams. Arguments in favor of higher order interactions. The American Naturalist, 121(6))8878811, 1983.
[6] Margaret Mayfield and Daniel B Stouffer. Higher-order interactions capture unexplained complexity in diverse communities. Nature ecology & evolution, 1(3):1-7, 2017.
10

---

[7] Shan Yu, Hongdian Yang, Hiroyuki Nakahara, Gustavo S Santos, Danko Nikolic, and Dietmar Plenz. Higher-order interactions characterized in cortical activity. Journal of neuroscience. 3 11(48):7554-11556 2011.
[8] Elad Schneidman, Michael Berry, Ronen Segev, and William Bialek. Weak pairwise correla- tions imply strongly correlated network states neural population. Nature, 440(7087):1007-
1012, 2006.
[9] Michael T Schaub, Austin R Benson, Paul Horn, Gabor Lippner, and Ali Jadbabaie. Random walks on simplicial complexes and the normalized hodge I-laplacian. SIAM Review, 62(2):353-
391, 2020.
[10] Timoteo Carletti, Federico Battiston, Giulia Cencetti, and Duccio Fanelli. Random walks on hypergraphs. Physical review E, 101(2)0223088 2020.
[1I] Unai Rldrrez- Federico Battiston, Guilherme Ferraz Arruda, Yamir Moreno, Matjaz Perc, and Vito Latora. Exolutionary dynamics of higher-order interactions in social networks. Nature Human Behaviour, 5(5):586-555, 2021.
[12] Ana Millan, Joaquin Torres, and Ginestra Bianconi. Synchronization in network geometries with finite spectral dimension. Physical Review E, 99(2):022307, 2019.
[13] Alexis Arnaudon, Robert L Peach, Giovanni Petri, and Paul Expert. Connecting hodge and sakaguchi-kuramon through mathematical framework for coupled oscillators on simplicial complexes. Communications Physics, 5(1):211, 2022.
[14] Charlotte Luff, Robert I. Peach, Emma-Jane Mallas, Edward Rhodes, Felix Laumann, Edward S Boyden, David J Sharp, Mauricio Barahona, and Nir Grossman. The neuron mixer and its impact human brain dynamics. bioRxiv, pages 2023-01, 2023..
[15] Marco Nurisso, Alexis Arnaudon, Maxime Lucas, Robert L. Peach, Paul Expert, Francesco Vaccarino, and Giovanni Petri. A unified framework for simplicial kuramoto models, 2023.
[16] John White, Eileen Southgate, J Nichol Thomson, Sydney Brenner, et al. The structure of the nervous system of the nematode caenorhabditis elegans. Philos Trans R Soc Lond B Biol Sci, 314(1 135)) 1986.
[17] Ronald Harry Atkin, Mathematical structure in human affairs. Heinemann Educational London, 1974.
[18] Jacopo Grilli, Gyorgy Barabas, Matthew J Michalska-Smitt, and Stefano Allesina. Higher-order interactions stabilize dynamics in competitive network models. Nature, 548(7666)2100-13,, 2017.
[19] Jean-Gabriel Young, Giovanni Petri, and Tiago P Peixoto. Hypergraph reconstruction from network data. Communications Physics, 4(1):1-11, 2021.
[20] Huan Wang, Chuang Ma, Han-Shuang Chen, Ying-Cheng Lai, and Hai-Feng Zhang. Full recon- struction of simplicial complexes from binary contagion and ising data. Nature communications, 13(1)33043, 2022.
[21] Andrea Santoro, Federico Battiston, Giovanni Petri, and Enrico Amico. Higher-order organiza- tion of multivariate time series. Nature Physics, pages 1-9, 2023.
[22]
H.O.Lancatter. The Chi-Squared Distribution. Wiley, 1969.
[23] Aleks Jakulin and Ivan Bratko. Quantifying and visualizing attribute interactions: An approach based on entropy. A1/30005. 04 2004.
[24]
Anthony J Bell. The co-information lattice. In Proceedings of the fifth international workshop on independent component analysis and blind signal separation: ICA, volume 2003, 2003.
[25] Fernando E Rosas, Pedro AM Mediano, Michael Gastpar, and Henrik J Jensen. Quantifying high-order interdependencies via multivariate extensions of the mutual information. Physical Review E, (0003)):33305,22011

---

[26] Bernd Streitberg. Lancaster interactions revisited. The Annals of Statistics, pages 1885, 1990.
[27] Edward H Ip, Yuchung J Wang, and Yeong-nan Yeh. Structural decompositions of multivariate distributions with applications in moment and cumulant. Journal of multivariate analysis, 19-1),119-134 2004.
[28] Arthur Gretton, Kenji Fukumizu, Choon Teo, Le Song, Bernhard Schölkopf, and Alex Smola. A kernel statistical test of independence. Advances in neural information processing systems, 20, 2007.
[29] Kacper Chwialkowski and Arthur Gretton. A kernel independence test for random processes. In International Conference on Machine Learning, pages 1422-1430. PMLR, 2014.
[30] Felix Laumann, Julius von Kügelgen, and Mauricio Barahona. Kernel two-sample and indepen- dence tests for nonstationary random processes. Engineering Proceedings, 5(1), 2021.
[31] Felix Laumann, Julius von Kugelgen, and Mauricio Barahona. Non-linear interlinkages and key objectives amongst the paris agreement and the sustainable development goals. arXiv preprint htt::2004.09//15 2020.
[32] Niklas Pfister, Peter Buhlmann, Bernhard Scholkopf, and Jonas Peters. Kernel-based tests for joint independence. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 80(1))5-31, 2018.
[33] Zhaolu Liu, Robert L. Peach, Felix Laumann, Sara Vallejo Mengod, and Mauricio Barahona. Kernel-based joint independence tests for multivariate stationary and nonstationary time-series. arXiv preprini arXiv:2305.08522 2023.
[34] Dino Sejdinovic, Arthur Gretton, and Wicher Bergsma. A kernel test for three-variable interac- tions. Advances in Neural Information Processing Systems, 26, 2013.
[35] Paul K Rubenstein, Kacper P Chwialkowski, and Arthur Gretton. A kemel test for three-variable interactions with random processes, arXiv preprint 1rX3v:/603.00922 2016
[36] Garrett Birkhoff. Lattice theory, volume 25, American Mathematical Soc., 1940.
[37] TP Speed. Cumulants and partition lattices 1. Australian Journal of Statistics, 25(2):378-388,
1983.
[38] Ginestra Bianconi. Higher-order networks. network.. Cambridge University Press, 2021.
[39] S-I Amari. Information geometry on hierarchy of probability distributions. IEEE transactions on information theory, 1715)117011171,, 2001.
[40] Mahito Sugiyama, Hiroyuki Nakahara, and Koji Tsuda. Tensor balancing on statistical manifold. In International Conference on Machine Learning, pages 2770-3279, PMLR, 2017.
[41] Edward H lp. Yuchung J Wang, and Yeong-nan Yeh. Some equivalence results concerning multiplicative lattice decompositions of multivariate densities. Journal of multivariate analysis, 84(2):403-409, 2003.
[42] Arthur Gretton, Karsten M Borgwardt, Malte J Rasch, Bernhard Scholkopf, and Alexander Smola. A kernel two-sample test. The Journal of Machine Learning Research, 13(1):723-773, 2012.
[43] George A Gratzer. Lattice theory: foundation, volume 2. Springer, 2011.
[44] Peter McCullagh. Tensor methods in statistics. Chapman and Hall/CRC, 2018.
[45] Eric Temple Bell. The iterated exponential integers. Annals of Mathematics, pages 539-557, 1938.
[46] Miklos Bona and Istvan Mezo. Real zeros and partitions without singleton blocks. European Journal of Combinatorics. 51:500-510, 2016.
12

---

[47] Lam Chi-Chin -h-nCC P Sadayappan, and Rephael Wenger. On optimizing a class of multi- dimensional loops with reduction for parallel execution. Parallel Processing Letters, 7(02):157- 168, 1997.
[48] Robert NC Pfeifer, Jutho Haegeman, and Frank Verstraete. Faster identification of optimal contraction sequences for tensor networks. Physical Review E, 90(3):033315, 2014.
[49] Frank Schindler and Adam Jermyn. Algorithms for tensor network contraction ordering. Machine Learning: Science and Technology, 1(3):035001, 2020.
[50] David C Van Essen, Stephen M Smith, Deanna M Barch, Timothy EJ Behrens, Essa Yacoub, Kamil Ugurbil, WU-Minn HCP Consortium, et al. The WU-Minn Human Connectome Project: An overview. Neuroimage, 80:62-79, 2013.
[51] Matthew Glasser, Stamatios Sotiropoulos, Anthony Wilson, Timothy Coalson, Bruce Fischl, Jesper I. Andersson, Junqian Xu, Saad Jbabdi, Matthew Webster, Jonathan R Polimeni, al. The minimal preprocessing pipelines for the Human Connectome Project. Neuroimage, 80:105-124, 2013.
[52] Alexander Schaefer, Ru Kong, Evan M Gordon, Timothy O Laumann, Xi-Nian Zuo. Avram J Holmes, Simon B Eickhoff, and BT Thomas Yco. Local-global parcellation of the human cerebral cortex from intrinsic functional connectivity MRI. Cerebral Cortex, 28(9):305-31144 2018.
[53] BT Thomas Yeo, Fenna M Krienen, Jorge Sepulcre, Mert R Sabuncu, Danial Lashkari, Marisa Hollinshead, Joshua L Roffman, Jordan W Smoller, Lilla Zollei, Jonathan R Polimeni, al. The organization of the human cerebral cortex estimated by intrinsic functional connectivity. Journal of Neurophysiology, 2011.
[54] Andrea Luppi, Pedro AM Mediano, Fernando Rosas, Negin Holland, Tim D Fryer, John O' O'Brien, James B Rowe, David K Menon, Daniel Bor, and Emmanuel A Stamatakis. A synergistic core for human brain evolution and cognition. Nature Neuroscience, (5(6)777 2022.
[55] Thomas F Varley, Maria Pope. Maria Grazia Puxeddu, Joshua Faskowitz, and Olaf Sporss, Partial entropy decomposition reveals higher-order structures in human brain activity. arXiv preprint arXiv:2301.05537, 2023.
[56] Paul L Williams and Randall D Beer. Nonnegative decomposition of multivariate information. arXiv preprint arXiv: 2010.
[57] Robin AA Ince. The partial entropy decomposition: Decomposing multivariate entropy and mutual information via pointwise common surprisal. arXiv preprint arXiv: 2017.
[58] Lennart Van Hirtum, Patrick De Causmaecker, Jens Goemaere, Tobias Kenter, Heinrich Riebler, Michael Lass, and Christian Plessl. A computation of d (9) using fpga supercomputing. arXiv preprint 20X.:2230440303,, 2023.
13

---

Appendix to "Interaction Measures, Partition Lattices and Kernel Tests for High-Order Interactions"
A Proofs
A.1 Proof of Proposition 1
First we prove the following:
Lemma 2. P I..d 1, for arbitrary i, then DiP 0
Proof. Without loss of generality let for expression partitions. For example, partition {(x², xx)) {x²), {X+)} can be written as 13/2/4 We denote the partition 1|23. as и. Notice that the factorisation only consists of an extraction of singleton, i.e. factorise one variable out of the joint distribution. The Möbius coefficients [4] for partial partition lattice above are.
simplify the notation, we adopt a shorthand convention
mla, 1) 1-1)--1(d-11
i
((11) 01-1
otherwise
The partial partition lattice is constructed
123 .
1|23.... ... d
2
...
12.....
1/2|3... 1/3/2....
2/3 ...d
1|2]3 ld
Given F1... = P.P2-... the partial partition lattice associated with Lancaster interaction collapses into the following:
1 1/23 d 2/1 ..d d
1 .....
3/122... 2/3
1/2/3...
1 13/2.....
1 ... |d
partition a is not refinement m, then it will be transformed to ar n, a partition one level down the lattice (as shown in the lattice above), e.g., let a be 13. d, then is be 2|1 ... d. When a a refinement of n, then
=
The transformed partitions will be cancelled out with the refinements of due to the fact that all partitions are counted exactly once and the partitions two consecutive levels have alternating signs. On the (d - 1) level of the lattice, exactly (n 1) partitions are not refinement of which are just enough to cancel out (n 1) number of ( at the bottom level.
When P II... Pr, where T. only contains non-singleton blocks, the Lancaster interaction
IP
II A
BETТ.
which is not necessarily zero unlike the Streitberg interaction. It will be zero, however, if any of the A is zero. If we factorise singleton of arbitrary block b, m, will be partition with one singleton block and hence A&P 0 Lemma above and does A P. Obviously if more singletons factorised out, A P will remain zero.
14

---

A.2 Proof of Proposition 2
Notice that this simplified Mixed Central Moment operator expression resembles the original Lan- caster interaction in Equation (2) using the unique ** notation' Notice that the resulting terms after expanding Equation 2 are distinctive with coefficient (-1) 1) | depending on the except we have multiple occurrences of Ô. This is due to the fact that PF; is the same as P.; For example, when d = 4, instead of direct product of the marginals P1 31 that accounts for Ó, there are also four other terms 1. P1 F1P51P33Ps P2133P 4. 1P 121PP that are also equal to Ô. Note that these terms
have different signs comparing with F2 P2IP3I because of P, and therefore result in (-1)d- 1) in the Möbius function.
A.3 Proof of Proposition
Computed immediately from the estimated norm of Mixed Central Moment operator, the Lancaster interaction estimator is
n
Malliss
II
(ab))
ΣΣЛ
(wr)) (<t))
ΣΣ ====li=
ab.
where Ø(.) = x()) - E[o(-)]-
A.4 Proof of Proposition
This can be completed by manual expansion and then equating the terms. Alternatively we can prove it more easily by the relationship of Lancaster interaction and Streitberg interaction Lemma] [
M₁ = E
II igb
E[o"]]
Kd
Σ ETE(ID)
(|TT)) 1)1(-1)(=11-
II M₁
bERT
Σ (|TAl 1)((-1)
E
Il iEb
(0 - E[0']]
{.EII(D)
bETT
A.5
Proof of Proposition S
Similar to the norm of the Mixed Central Moment operator, the norm of Mixed Cumulant operator can be obtained by simply squaring. Note that this can be computed directly using the product lattice discussed in Section 4.2
[alll
(|TTLl - 1)1(-1)(*~11-
EII
お (xt.),
1)((-1)(=211-
3
(afk)
n M
{b,} ==1
[b/1 }
i=
Σ
(|T
1)!(1=') ()((-))(F=)-))) (1+
1
ΣII
3d (ot.). Σ
II
35 (x
n. m7
(b.) =11
(0))
Σ
(|T - 1)((=(]] - 1)((-1) 1)((-1)=01++)))
ΣΣ
K
n n
(b.) 10))= 1=1
15

---

B Product Lattice
The partition lattice is reduced to upper semilattice after centring. When d = 4, the reduced lattice only contains 12234, S12 14 and P14 which are the 4 partitions with no singletons. Given the two arbitrary lattices L and L' the product lattice formulated using the Cartesian product of S x S" with product order can be defined. Given two elements (a, o'), (II, x²) E x L', (a, & (л,, x²) if and only if a S and o' I '². In our case specifically, L and L' are the identical reduced partition lattices of variables after centring, and the elements in the product lattice are analogous to the inner products when computing the kernel-based estimator. Below show the product lattice Therefore, it can be utilised to directly compute the coefficients of the estimator in Proposition 5 Note that the order within each element not particularly informative in this case as the kernel matrices are symmetric.
4.
Product lattice 4
X
Figure 4: Product lattice for d 4. The lines and squares represent the non-singleton factorisations with two blocks and one block. Dark green and dark red highlight the elements in C, whilst light green and light red are used for the elements in L'.
The resulting estimator using centred kernels when d = 4 is
|14
K
K
Kt
2
2
2
K
n
Kjk
Klykg
Σ
K Σ
ijkt
jjkl
GR
K
K
Σ
Σ
n
Kkk
Notice that a few terms are summed up twice. This is due to the symmetry in the product lattice we described above.
C Composite Null Hypothesis
The permutation test strategy we propose for Streitberg interaction and Lancaster interaction is carried out via multiple sub-hypotheses, and the rejection of the interaction is confirmed once all the sub-hypotheses are rejected. Instead of having to consider all possible factorisations (Ba) for sub-hypotheses, by considering the lattice structure we only need to perform (2 1) sub-tests of factorisations corresponding to the 2nd level of the partition lattice, since all other sub-tests are consequences of these.
In the case of rejecting the Streitberg interaction, won't necessarily be the case that all sub- hypotheses are rejected, allowing us to derive a more detailed understanding of how the joint distribution is factorised. For example, in the 7 sub-hypotheses for the Streitberg interaction when d = 4, there is a clear difference between the number of rejections when the ground truth factorisations are P, 234 or PP P2I234. The former results rejections while the latter results in only rejections because IP PPP is the mutual refinement of IP12PPa4 34, IP, P224 and PPP1344 Hence, even though the vanishing conditions of Streitherg interaction are not 'if and only if", we can still narrow down the range of the possible ground truth factorisations using the rejected sub-hypotheses.
16

---

The operations can be easily explained using the lattice. Before performing any sub-test, II(D) the space of all possible ground truth factorisation. Whenever we reject a sub-hypothesis, the refinements of the associated partition и gets eliminated. For example, if we reject IP] P234. then P2 and P/ In other words, partitions in the interval lattice [Ū, r are eliminated. The remaining partitions until there are no further rejections are the only possible choices and also form a lattice. Again we illustrate this idea using the partition lattice variables (see Figure 5
also automatically reject P1 3422334 P1 P1P31P24.
Figure 5: Remaining factorisations with partial rejection of sut-heyot for d = 4. The greyed partitions represent the eliminated factorisations due to the rejections in the sub-tests. A dashed line between two partitions indicates the refinement ordering (only shown for those not rejected). automatically
(a) In the case that 4 of the sub-hypotheses are rejected, all of their refinements climinated. The remaining elements form an interval lattice and is isomorphic to the partition lattice of variables. (b) If 6 of the sub-hypotheses can be rejected, the ground truth can only be either the factorisation remaining or non-factorisable. Similarly this is isomorphic to the partition lattice of two variables.
Remark: The simplification of the operators (computation of the test statistics) can be reflected in collapsed partition lattice without singletons. However, one should not confuse this with the full lattice that is associated with the test procedures, i.c. one should always consider the full partition lattice to figure out the possible configurations of the ground truth factorisation based what we outlined above. The two lattices have the same form however they have completely different implications. The lattice corresponding to test procedures only be reduced result of the rejections. For example, sub-tests are needed instead just 3 tests that involve the factorisations without singletons when 4.
D
Möbius Inversion
Instead of deriving the Möbius function as the inverse of the Zeta function, it can also be computed using the expression below (37]
a= π o <<T. otherwise
((a,z) = M(s,o)
0
We can check the validity of this formula by computing the Möbius on the trivial partition lattice associated with joint independence, which contains two elements, Ô and i,
(1i, 1)
u 1) (mi, -1,
hence
A/p p(i, ((,,))PI.... I..d + (0, 1)
d
P i
P,
1=1
17

---

E Links to Simplicial Complexes
k-simplex is set of k+1 vertices 0 [po 1Poo,..,,Ph pk], i.e., 1-simplex is line, 2-simplex is triangle, 3-simplex is a tetrahedron, etc. simplicial complex @ is a collection of simplices that satisty two conditions: (i) if θ, then all the sub-simplices contained in θ; and (ii) the non-empty intersection of two simplices 0,, 02 E θ is a sub-simplex of both 0 and 0z. The first condition makes the inclusion ordering in the subset lattice natural for simplicial complexes.
0 built from subsets of 0 also
Below we have the Zeta matrix in Figure 6 and Möbius matrix in Figure 7 illustrate that the boundary operators (with no directionality) of (d-l)-simplex can be founded as block matrices within them. The subset lattice doesn't include the non-simplicial terms and we have shown numerically that these terms are essentially in the synthetic experiment when P1234 P12IP34
13/34
1/1//4
1254
125
134
284
13]24
14/23
34
E-F
Non-simplicial
Figure 6: Zeta matrix for the partition lattice of 4 variables. Singletons in the partial factorisations are omitted for simplicity. The matrix can be decomposed into small blocks which correspond the boundary operators in 3-simplex. Here the two red blocks represent the Edge-to-Face operator and Terraho-Ter operator. The blue block corresponds the non-simplicials, partitions with no singletons.
1284
128
124
134
234
12/34
18/24
14/28
12
18 14 23 24
24 1/2/3/4
1234
F-T
124
131
1234
1524
-1
1/23/4
.8
E-F
Non-simplicial
Figure 7: Möbius matrix for the partition lattice of 4 variables. Comparing with the Zeta matrix above we see that the block structure preserved. Each block can be ob- tained using the Möbius inversion e.g. {12]34, 13/24, 14/23, 12/3/4, 13/2|4, 14/2/3, 23) 1|4, 24|1|]334|||2] forms a partially ordered set. The values are always -1 since the partially order set only contains the elements from two neighbouring levels.
the incidence matrix of its partially ordered set,
18

---

F
Streitberg-> Interplay of Test Power and Sample Size
To explore the effect of the different sample size on test power, we have performed four experiments, each with three interaction strengths, to show how decrease in the sample size affects the null rejection rate (Fig In all cases, the Lancaster and Streitberg tests are able to accurately detect higher order interactions with low 50- 100 samples (depending the interaction proportion//srennthh), whilst dHSIC requires substantially more samples. This provides further experimental evidence that the proposed Streitberg interaction test has superior sensitivity, compared dHSIC, in detecting higher order interactions.
4-var independence dHSIC lencasior Stroiberg
4- 4-way interaction
5-var independence
5-way interaction
aPo Stre trarg
modfed
Strviberg
rolectoon rate
1.00
1.00
1.00
1.00
100%
0.75
0.75
0.75
0.75
0.00
Interaction proportion
10
102
10'
10/
10'
10'
30
10'
10*
null relection rate
1.00
1.00
1.00
1.00
75%
0.75
0.75
0.75
0.75
0.50
0.50
0.50
0.50
10
102
10
10a
10
10a
10a
10
102
107
mull rejection rate
1.00
1.00
1.00
1.00
50%
0.75
0.75
0.75
0.75
0.50
0.50
0.50
0.50
0.25
0.25
10
10'
10
10'
10
10"
10*
10
10'
10a
sample Size
sample size
sample 528
stmple
Figure 8: Decay of accuracy rates as the sample size is decreased in the XOR example. Inde- pendence tests (first and third columns); interaction tests (second and fourth columns). In all cases, dHSIC degrades substantially faster than Lancaster and Streitberg.
G Computational Complexity and Practical Limitations
G.1
Computational Complexity
The original time complexity for computing the Streitberg interaction estimator is O(3rn2d where Ba is the Bell number, which represents the number of partitions in a set of cardinality d. Hence B3 the number of terms in the Streitberg interaction estimator. For fixed number of samples n, our lattice formulation allows to reduce the number of terms in the mixed cumulant operator from Ba to Fa. Table Table] iilustrates this reduction in the number of terms before and after eliminating the partitions with singletons. After this reduction, the complexity becomes (()),2dd)
Table Number terms in the test statistics before/after eliminating the partitions with singletons.
No. of variables
4
5
6
7
8
15->(4)
52->(11)
201->(41)
877->112)
4140=>7755) 248->(1)
Laneaster=>(eentred)
12=>(1)
27=(())
58->(1)
121=(1)
For fixed d, we can further optimise the time complexity by optimal contraction ordering, using the fact that two index sets are disjoint, such that (n) becomes min,ii((() This reduction for different d can be found in T
19

---

Table 2: Time complexity for d-order interaction estimators.
d-order
dHSIC
Lancastrr=Xopimmiee))
Streitberg-Xoptimissdd
2-way
(n²))
O(n")
O(n²)
3-way
(nn))
O(n")
O(n))
4- 4-way
On
O(n")=
5- 5-way
O(n²)
=0())=((m))
)(n²)=
Similarly, the time complexity for computing the d-order Lancaster interaction estimator naively is O(2 ?n) By centring, i.c. eliminating the partitions with singletons, the only element left is 1. Therefore the Lancaster interaction estimator can be computed in O(d"))
G.2 Practical Considerations
The problem of detecting d-order interactions among a group of e variables combinatorial, as it entails checking groups of d variables chosen from the variables. Clearly, such combinatorial problems become infcasible as e and/or d become large. Recent work on real data from various application areas has focused shown that interactions with d = 3, 4 already make a significant difference to the analysis of network structure and network dynamics 2 21 (55, highlighting the potential benefits to be gained from tests that detect high order interactions.
revealing any interactions beyond pairwise, i.c., d > 2. It has been
Many well-recognised and widely used theoretical approaches only have closed forms for d 3 [56 57 and cannot be generalised arbitrary d, since the number terms those approaches related to the Dedekind number which becomes rapidly intractable (58) In contrast, here we show that the Streitherg interaction can be explicitly defined for any d, and we have devised theoretical and
computational strategies to reduce its computational cost via the lattice theory formulation.
In the case where high order interactions for range of d are of interest, possible to discount the computation of certain high order interactions when some lower-order interactions are present. If all the lower order interactions are absent, then testing the d-order interaction is the same as testing the joint independence for d variables. For example when = 3 and all pairwise interactions are absent, then testing the Streitberg interaction reduces testing joint independence of variables. More generally, if we test the interactions bottom-up, i.e., recursively from the lower orders upwards, the expression for Streitherg becomes simpler whenever there is a lower-order independence. The simplification is possible because the Streitberg interaction can be rewritten as the sum of the differences between a factorisation and the product of the marginals due to the fact that the sum of Möbius coefficients is zero. Hence the presence of lower order independence allows to simplify the d-order interaction formula.
Alternatively, there are also cost- nedtreeu simplifications if we do the tests top-down (i.e. starting from order d downwards), although the problem still remains combinatorial. there are partial rejections for some tests involved in the d-order Streitberg test, we can narrow down the possible choices of the factorisation (as discussed in Appendix C). This allows us to eliminate the lower order factorisations that fall in the lattice branches of the rejected second level factorisations, thus reducing the total tests needed for the factorisation of the d-variables. Further simplifications are achieved by accounting for overlaps of the lattice branches.
All these observations can be taken into account in the construction of hypergraph representations, as further discussed below in Appendix H
H Neuroimaging Data
H.1 Preprocessing
Preprocessing of neuroimaging data was performed following Luppi et al. (54]. We summarise the key steps here, and refer interested readers to the original paper for further details.
We started from the 'minimally preprocessed" release of the fMRI data from the Human Connectome Project [50 51]. This data preprocessed with bias field correction, functional realignment, motion
20

---

correction and spatial normalization to Montreal Neurological Institute (MNI-152) standard space with mm of isotropic resampling resolution. We then removed the first ten points in the time series to avoid transient effects introduced by the scanner. Finally, we further denoised the data with the anatomical CompCor method, which involves regressing out potential noise confounds (specifically, five principal components of white matter activity, five principal components of cerebrospinal fluid activity, and 12 motion parameters including head translation, rotation, and their temporal derivatives). All preprocessing steps were performed using the CONN toolbox (https: www org. ptrg/projects/. cetojjecs/coon////// ), version 17158. The resulting volume was parcellated according to the Schaefer- 100 atlas [52] by spatially averaging across all voxels in the same region for each timestep.
Table 3: Time taken for experiments in Figure 3
SOM
VIS
SAL
DAN
DMN
FPN
LIM
Random
2-way
1s
2s
1s
2s
2s
2s
1s
12s
3-way
12s
18s
7
13s
8s
8s
1s
13s
4-way
5m36s
5m
2m44s
2m31s
2m10s
1s
2m41s
5-way
2h18m24s
1h30m9s
1h3m
54m40s
1h 1m44
27m16s
2s
29m
We report the computational time for all fMRI experiments in Table All experiments carried
a 2015 iMac with GHz Quad-Core Intel Core i7 processor and 32 GB 1867 MHz DDR3 memory.
H.2 Analysis of Simple Hypergraphs
As an illustration of future lines of work, we hierarchically constructed hypergraphs from the identified d-order interactions. Specifically, individual regions nodes and the high-order interactions used to define hyperedges, incrementally including the interactions of increasing order. The resulting hypergraphs integrate information across orders of interaction, and can be analysed by computing structural properties of hypergraphs. As an example, we show the degree assortativity of the different hypergraphs in Fig 9 We find that almost all RSNs display negative degree assortativity (but much less so than random), and only FPN shows positive degree assortativity for higher order hypergraphs. Future analyses of these hypergraphs will study the links between high-order interactions in brain activity and different functional areas.
SON
2.8
DAN
0.2
auartulvkyy
6.
Ran-ie
-0.2
-0.6
0.0
-0.0
2-way
3-way
4x0ay
5-wayy
2
3+2 H00
4+3+2wMy
5444342-way
Figure 9: Hypergraph analysis of the neuroimaging data. Percentage of high-order interactions in RSNs (left) and degree assortativity of hypergraphs constructed from the high-order interactions detected (right).
H.3
Alternative Hypergraphs: Emergent and Redundant High- Interactions
To interpret the presence of lower order interactions in larger cliques, one can alternatively build hypergraph as follows: i) If a d-order hyperedge is detected and none of the lower order (d-1)-order hyperedges are present then we say that this d-order hyperedge reflects a purely synergistic (or emergent) interaction between the d variables; ii) if the d-order and all the lower order interactions are present, then we say that the d-order interaction is purely redundant (and corresponds a simplicial complex construction); iii) if some, but not all, lower order interactions are present (e.g., for 4 variables, only the 4-way interaction and one 3-way interaction are present), then both synergy and redundancy present in this group of variables. Such construction could offer alternative.
21

---

statitically-motiivee motivated approach to computing synergy and redundancy, an important current topic of research in computational neuroscience, and its representation through hypergraphs.
I Code
This code for performing the interaction tests and the synthetic experiments is provided in this anonymous Github repository ltttps://github cow/barahonaaaaa..aaaaaaaaaaaaaa - wwwpeahhhaarrrraceeeeeaaaae www.erre-eerrreeeeeeeeeee group/ streitberg iat.beergginteracciin......
22


Reviews:
[{"review_id": "iX2dNRWzf8", "reviewer": "Reviewer_SriA", "summary": "The paper is a nice introduction to measuring high order interactions between groups of random variables and an experimental demonstration of measuring these interactions for up to 4 variables inclusively. The main focus is on the Lancaster and Streitberg interactions tests via kernel methods such as the Hilbert Schmidt Information Criterion. After presenting a mix of prior and, possibly, novel results some experimental evidence is provided that the approach works.\n", "questions": "1. What is the single technical novel result in this paper?\n2. Experiments demonstrating the effect of the interplay between the accuracy and computational complexity. Computational complexity numbers reported for all runs of the fMRI experiments.\n3. Error bars in Figure 3 are missing. Are results so stable?\n(the authors have responded to my questions, which is now reflected in my raised score)", "limitations": "nothing to report", "rating": 5, "confidence": 2, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "1. The topic of capturing interactions beyond binary is interesting and largely avoided for lack of adequate tools. \n2. The paper is generally well written and the didactic value is great. It can serve as an accessible introduction into the topic.\n", "weaknesses": "1. The paper's contribution is unclear. Some of the claims from the abstract are already known previously published work, such as the link to the lattice theory (Streitberg) and the work with d=4 [33, 34].\n2. The paper claims efficiency of the described permutation procedure but besides high level considerations of its asymptotic behavior no experimental characterization is provided. Especially the interplay of efficiency and accuracy. When the procedure is supposed to be faster for a smaller number of samples how fast does the accuracy drop?\n3. I found the presentation unclear. On one hand extreme novelty and generality with respect to d is claimed. On the other, experiments are only shown for d up to 4 (previously available results).\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their comments. \n\n- **Novelty**: Although, as the reviewer points out, the mathematical link between lattice theory and the $d$-order probabilistic Streitberg measure was indicated in [25] and briefly mentioned in [33], in neither reference were interaction tests implemented owing to their computational complexity. (Reference [34] does not make any references to lattice theory.) Regarding novelty, the main theoretical contribution here is to leverage lattice theory to formulate and optimise $d$-order Lancaster and Streitberg interaction tests that are computable on real data, thus enabling practical applications. More specifically:\nFirst, we formulate the differences between joint independence, Lancaster interaction and Streitberg interactions in terms of the partition lattice and its sublattices. \nSecond, we define the $d$-order Streitberg and $d$-order Lancaster interaction tests (with improved understanding of the vanishing conditions on the Lancaster interaction) in the kernel setting, and optimise their computational efficiency by eliminating partitions with singletons in the respective lattices.\nThird, we propose a recursive computation strategy such that one would only need to compute the inner product if the join of the two partitions is $\\hat{1}$.\nFourth, after the optimisation using centring, the Hilbert-Schmidt norm of each high-order measure can be directly obtained from the product lattice (explained in detail in section B of the SI).\nFinally, we show that using the interval lattice we can derive a generalised interaction measure whose test formulation can be directly obtained from the 2nd level of the lattice which captures the corresponding factorisations in the sub-hypotheses.\n\n- **Computational efficiency**: We apologise for our lack of clarity. Our comment about efficiency relates to the fact that instead of having to consider all possible factorisations ($B_d$, Bell number of $d$ variables), we only need to perform $(2^{d-1}-1)$ sub-tests of factorisations corresponding to the 2nd level of the partition lattice, since all other sub-tests are consequences of these. This is a substantial reduction in computation that follows from the lattice formulation. We have corrected our wording in the main text to make this clearer.\n\n- **Interplay of efficiency and accuracy**: We thank the reviewer for this remark. To explore this, we have performed 4 different experiments each with 3 different interaction strengths to show how decreasing the sample size affects the null rejection rate (Fig 1, attached). In all cases, as the sample size decreases, the null rejection rates of dHSIC decay much faster than both Lancaster and Streitberg.\n\n- **Higher order $d>4$ and practical applications**: To highlight the generality of the proposed statistical tests to any order, we have performed $d=5$ numerical experiments for both the synthetic examples (Fig 2, attached) and the real-world data (Fig 3, attached). The decay of null rejection rate with decreasing sample size for $d=5$ is included in Fig. 1.\n\n- **Our empirical results provide backing for  practical applications**: We have shown that: (i) Lancaster is better than dHSIC for testing joint and marginal independence; and (ii) Streitberg is more data-efficient compared to dHSIC for detecting high-order interactions. Moreover, we show that these tests are practical to compute on standard computers for iid data, and we have now added further examples of applications to stationary time-series data (Fig 4, attached).\n\n- **Computational complexity**: The computational complexity of dHSIC, Lancaster and Streitberg tests for $d=2,3,4,5$ are shown in Table 1 in the main response. For $d=4,5$, we show both the full complexity of the na\\\"ive enumeration compared to the substantially reduced complexity after optimisation using the lattice formulation. We also show the computational time for all fMRI experiments in Table 2. All experiments carried out on a 2015 iMac with 4 GHz Quad-Core Intel Core i7 processor and 32 GB 1867 MHz DDR3 memory.\n\n- **Error bars**: Thank you for spotting this error. Error bars in the bar chart now added to Figure 3. The results are stable, and this will be changed in the revision."}, {"review_id": "vOHXtypVGc", "reviewer": "Reviewer_KEns", "summary": "The authors conduct a very systematic study of tests that measure couplings between variables where these couplings include higher-order interactions. Authors uncover connections with lattice theory that then leverage them to formulate these tests more efficiently and in a more interpretable way. The authors demonstrate empirically that one of these tests, the Lancaster test, is better than dHSIC to test for joint and marginal independence and that the Streitberg test is overall superior to the other tests when it comes to detecting all the factorizations of the joint distribution. Finally, the authors validate their results numerically with synthetic data and brain data.", "questions": "No questions for the authors.", "limitations": "The authors have addressed adequately all the limitations and negative impact of their work.", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- **originality**: The authors provide novel theoretical connections with lattice theory, which help derive the interaction measures and their statistical tests. These contributions are novel to the best of my knowledge and are noteworthy.\n\n- **quality**: The theoretical derivations seem sound and there are no obvious errors I could identify. The theoretical results are validated with a sufficient amount of simulated and real-world data.\n\n- **clarity**: the paper is very clearly written. \n\n", "weaknesses": "- **significance**: the significance of the approach is hindered by the fact that it is, as claimed by the authors, \"computationally expensive\". Additionally, the authors themselves also point out that the theoretical results rely on the assumption of iid data, which limits the applicability of the method to real data. ", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their review and comments. The theoretical results in our paper leveraging lattice theory enable us to define and implement kernel-based $d$-order Lancaster and Streitberg interaction tests, whose computational complexity was previously prohibitive, thus opening up their practical use for real-world applications.\n\nThe original time complexity for computing the Streitberg interaction estimator is $\\mathcal{O}(B_d^2 n^{2d})$ where $B_d$ is the Bell number which represent the number of partitions in a set of cardinality $d$ and, consequently,  $B_d^2$ represent the number of terms in the Streitberg interaction estimator. For fixed $n$, our lattice formulation allows us to reduce the number of terms in the mixed cumulant operator from $B_d$ to $F_d$ where $F_d$ is the number of partitions without singletons, equivalent to centring. In Table 3 in the main response, we illustrate the reduction in the number of terms before and after eliminating the partitions with singletons. After this reduction, the complexity becomes $\\mathcal{O}(F_d^2 n^{2d})$. For fixed $d$, we can further optimise the time complexity by optimal contraction ordering, using the fact that two index sets are disjoint such that $\\mathcal{O}(n^{2d})$ becomes $\\mathcal{O}(n^{min(|\\pi_s|, |\\pi_s'|) + 1})$. Examples of this reduction for fixed $d$ can be found in Table 1 in the main response. \n\nSimilarly, the time complexity for computing the $d$-order Lancaster interaction estimator na\\\"ively is $\\mathcal{O}(2^{d+1} n^{2d})$. By centring, i.e. eliminating the partitions with singletons, the only element left is $\\hat{1}$. Therefore the Lancaster interaction estimator can be computed in $\\mathcal{O}(dn^2)$. \n\nAdditional examples with stationary time series data: \nFollowing the reviewer's comment on $iid$ data, we show that our $d$-order interaction tests can be generalised to time-series data too. The tests introduced here can be used in conjunction with a permutation procedure that approximates the null distribution by shifting through the time observations following a recent paper [1]. We are thus able to test high-order interactions between $d$ stationary random processes.\nWe have implemented this strategy and show some results from synthetic data in Figure 4 in the attached pdf. Again, we see that the Lancaster test is able to detect factorisations that corresponds to the partitions with at least one singleton, and the Streitberg test is able to detect any factorisations (dHSIC failed in both a and b). This means that both tests have controlled type I error. In Figure 4c, we also see that the Streitberg and Lancaster tests reach full power with fewer time observations compared with dHSIC when testing joint independence. Similarly when testing high-order interactions in Figure 4d, Streitberg is more data-efficient compared with modified dHSIC.\n\n[1] Liu, Zhaolu, et al. \"Kernel-based Joint Independence Tests for Multivariate Stationary and Nonstationary Time-Series.\" arXiv preprint arXiv:2305.08529 (2023)."}, {"review_id": "to1a4Tsoll", "reviewer": "Reviewer_vitW", "summary": "The paper describes a kernel test for estimating d-th order interaction, and computing its significance using permutation test. The authors propose two tests based on Lancaster interaction and Sreitberg interaction. The authors demonstrate the effectiveness of these measures on simulated data with higher order interaction and on real data from fMRI measurements.", "questions": "- equation 4, do we need the inequality under the sum given the zeta function?\n- equation 5, what is phi_iz", "limitations": "The authors address the limitation of the work in section 7, e.g., around computational complexity.", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 4, "strengths": "The paper addresses an interesting and relevant problem since estimating higher order interaction provides invaluable insight into the structure and relations among multivariate systems. The paper is very well written; the authors motivates and describes the problem well, and discusses the proposed estimator in detail.", "weaknesses": "The paper only discusses the performance of the proposed measure against a kernel based test of independence rather than existing tests on Lancaster and Streitberg interaction.\n\nIn the experimental section it is shown the neuroimage dataset shows higher order interaction within a region than between region. It will e great to get some more insight on, e.g., the implication of higher 4-way interaction that 3-way in SOM and higher 3-way interaction than 4-way interaction in DMN. Similarly, the lack on higher order interactions in LLM despite 2-way interactions being prevalent.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their comments and positive response to our paper. \n\nRegarding the suggestion to compare against non-kernel implementations of Lancaster and Streitberg tests, we are not aware of any such alternative approaches. However, we would be open to exploring alternatives for comparing and validating our implemented tests if the reviewer has suggestions.\n\nAlthough our focus here was on the theory and the neuroimaging data was presented very concisely, we were also intrigued by the analysis of the this data set and its interpretations. If we consider the SOM or VIS, our results show a larger percentage of 2-way interactions compared to regions such as the FPN or DAN, suggesting that there is increased redundancy in the system relative to these other regions. The SOM and VIS are structurally coupled, modular sensorimotor processing regions, which could benefit from information redundancy to increase robustness. These findings are in line with recent work reporting higher redundancy in sensorimotor regions [1], although further investigation is needed to draw firm conclusions about brain function from these analyses.\n\nThe potential for future investigation of the high-order features in neural data is large but this fell beyond the limits of the current paper. As an additional investigation, we hierarchically constructed hypergraphs from the identified interactions at increasing order $d$, by defining the high-order interactions as hyperedges and individual regions as nodes, incrementally including the interactions at a higher order. We then analysed some structural properties of these hypergraphs, and report the degree assortativity of each hypergraph (Fig 3, attached). We found that most regions display a negative degree assortativity that became less negative with increasing order interactions.\n\nAdditionally, we noticed that the degree assortativities when only considering pairwise connections d=2 are quite similar across regions. Standard deviations of degree assortativities across RSNs are shown in the table below.\n| Degree assortativity | Standard deviation |\n|----------------------|--------------------|\n| 2-way                | 0.044              |\n| 3+2-way              | 0.161              |\n| 4+3+2-way            | 0.104              |\n| 5+4+3+2-way          | 0.059              |\n\nHowever, when including high-order interactions in the hypergraphs, we find that the degree assortativities of each region diverge from each other. We interpret this as tentative evidence that the structure of macroscopic brain organisation may differ substantially when taking into account interactions beyond pairs, highlighting the importance of methods like ours to reveal new insights from brain data when taking into account high-order interactions that could be of importance for functional processing.\n\nWe also thank the reviewer for spotting the type error in Equation 4. You are right that the inequality has been encoded in the zeta function already. We have corrected this in the revised manuscript.\n\nRegarding your question on Equation 5, we are assuming that by $phi\\_iz$ you are referring to $\\phi^i$? $\\phi^i$ is the feature map of variable $X^i$ and $k^i(x, x') = \\langle \\phi^i (x), \\phi^i (x') \\rangle$. We apologise for not making this clear in the text and have stated this more clearly in Definition 1 and Definition 2 in the revised manuscript.  \n\n[1]  Luppi, A.I., Mediano, P.A.M., Rosas, F.E. et al. A synergistic core for human brain evolution and cognition. Nat Neurosci 25, 771–782 (2022)."}, {"review_id": "7ewGFtfwZP", "reviewer": "Reviewer_6qiu", "summary": "The authors introduce a hierarchy of $d$-order interaction measures by introducing a family of tests based on factorizations of the joint probability distribution that generalize to any order $d$ and define non-parametric tests to establish the statistical significance of these d-order interactions. They link their approach to lattice theory which they use to reduce the computational complexity of our d-order interaction tests. They validate their method on a synthetic data set and through an application to a neuroimaging dataset.", "questions": "Can the authors comment more on the current practical limitations?\nCan the authors better flesh out the consequences and interpretations of their results for the HCP data, as a means of demonstrating the utility and motivation for their approach?", "limitations": "The authors could explicitly add a \"Limitations\" subsection.", "rating": 6, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The paper is very well written, covers an exciting topic and seems to introduce very original results that are of broad significance. The tests proposed and their relationship to lattice theory are fascinating and of some practical utility, in my opinion (although I must admit this work is outside my expertise).", "weaknesses": "As pointed out by the authors, for the time complexity of computing the Streitberg interaction estimator, the number of terms remains combinatorial. Yet clearly, as the authors show, it is possible to apply it to the HCP data. It would be useful to know the current practical limitations in terms of data size for the method proposed here.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their comments and useful suggestions. We appreciate that one must consider both the overall complexity and the practical implementation of our tests. Therefore, we have provided a more in-depth analysis of both. \n\nFirst, in Table 2 in the main response we provide the overall time taken to compute the full experiments at each order for each brain network. Whilst the overall time to compute the tests increases with interaction order, it was still entirely feasible to run this on standard computer hardware, particularly given that 3-, 4-, and 5-way interactions will be sufficient for most datasets.\n\nSecond, datasets often come with different sizes and so we further considered how sample size would affect our ability to identify high-order interactions. In the attached PDF, we Figure 1 shows how the rejection rate of our tests reduces as we decrease the sample size. The Lancaster and Streitberg tests are able to accurately detect higher order interactions with as low as 50-100 samples (depending on the interaction proportion/strength), whilst dHSIC requires substantially more samples. \n\nFinally, we report the big-O notation complexity in Table 1 in the main response for each interaction order. We show that our theoretical improvements, which follow from our mathematical links with lattice theory, reduce the time complexity significantly for $d>3$ interactions in both the Lancaster and Streitberg interaction estimators.\n\nRegarding the HCP data set, we believe there are some interesting interpretations and important implications. If we consider the SOM or VIS, our results show a larger percentage of 2-way interactions compared to regions such as the FPN or DAN, suggesting that there is increased redundancy in the system relative to these other regions. The SOM and VIS are structurally coupled, modular sensorimotor processing regions, that benefit from information redundancy to increase robustness. These findings are in line with recent work reporting higher redundancy in sensorimotor regions [1], although further investigation is needed to draw firm conclusions about brain function from these analyses.\n\nFor further investigation, we hierarchically constructed hypergraphs from the identified interactions at different orders and then analysed their structure. By computing the degree assortativity of each hypergraph we found that most regions display a negative degree assortativity that became less negative with increasing order interactions. Additionally, we noticed that the degree assortativities when only considering pairwise connections d=2 are quite similar across regions. However, when including higher-order interactions in the hypergraphs, we find that the degree assortativities of each region diverge from each other. This has important implications on how the RSNs are defined if one considers high-order interactions.\n\n[1]  Luppi, A.I., Mediano, P.A.M., Rosas, F.E. et al. A synergistic core for human brain evolution and cognition. Nat Neurosci 25, 771–782 (2022). "}]

Author rebuttal:
We thank the reviewers for their clear reviews and thoughtful questions. In addition to the specific answers in the responses to each reviewer, we would like to briefly address here three overarching themes that have appeared in the reviews: the computational complexity of our method; implications of the analysis of real-world neuroimaging data; and the novel contributions of this manuscript above and beyond prior work. We address each of these in turn:

- **Computational complexity**. As previously stated, our method scales combinatorially in $d$, as would be expected for a method that evaluates $d$-order interactions. Nonetheless, for finite $d$ values relevant in practical applications our method achieves substantial and measurable improvements in complexity (see Table 1 below).
This was achieved by three technical contributions: 
    1. elimination of partitions with singletons via kernel centring;
    2. formulation of an optimal contraction ordering; and 
    3. a recursive computational strategy. 

    These reductions follow from our mathematical formulation and links to lattice theory. Moreover, we run various experiments and highlight the practicality of implementing these tests with different sample sizes and also with different data types (i.e., stationary time series) which were previously unattainable without our theoretical improvements. Some of this work was included in the SI, and additional work is explained in the individual responses below. In the revised version of the paper, we will include further discussion of prior work and add experimental results to show explicitly that the ways to optimise using lattice theory are novel leading to substantial reductions in computation.
- **Implications of the results in real-world data**. The key take-away message from our experiments --- beyond the fact that the proposed quantities are tractable --- is that brain activity has substantial and heterogeneous high-order structure that can be revealed using our method. This high-order structure is non-trivial, in the sense that it cannot be predicted from pairwise interactions, and calls for new hypergraph analyses to understand brain function across scales. In the revised version of the paper we will report new experiments emphasising the importance of this high-order structure and discuss links with supporting neuroscientific evidence.

- **Novelty**. The paper's main novel theoretical contribution is the use of lattice theory to formulate and optimise $d$-order Lancaster and Streitberg interaction tests. The novel contributions include:
    1. We use the partition lattice and its sublattices to make a qualitative comparison between the information encoded in joint independence, Lancaster interaction and Streitberg interactions.
    2. We formalise the $d$-order Streitberg and the $d$-order Lancaster in the kernel test setting and optimise their computational efficiency by eliminating the partitions with singletons in the respective lattices.
    3. We propose a recursive computation strategy such that one would only need to compute the inner product if the join of the two partitions is $\hat{1}$.
    4. After the optimisation using centring, the Hilbert-Schmidt norm of each high-order measure can be directly obtained from the product lattice (explained in detail in section B of the SI).
    5. We show that using the interval lattice we can derive a generalised interaction measure. The respective test formulation is also related to lattice since the corresponding factorisations in the sub-hypotheses can be directly obtained from the 2nd level of the lattice.

**Table 1**. Time complexity for $d$-order interaction estimators.
| $d$-order | dHSIC              | Lancaster$\Rightarrow$(optimised)                      | Streitberg$\Rightarrow$(optimised)                     |
|-----------|--------------------|--------------------------------------------------------|--------------------------------------------------------|
| 2-way     | $\mathcal{O}(n^2)$ | $\mathcal{O}(n^2)$                                     | $\mathcal{O}(n^2)$                                     |
| 3-way     | $\mathcal{O}(n^2)$ | $\mathcal{O}(n^2)$                                     | $\mathcal{O}(n^2)$                                     |
| 4-way     | $\mathcal{O}(n^2)$ | $\mathcal{O}(n^8)$$\Rightarrow${$\mathcal{O}(n^2)$}    | $\mathcal{O}(n^8)$$\Rightarrow${$\mathcal{O}(n^3)$}    |
| 5-way     | $\mathcal{O}(n^2)$ | $\mathcal{O}(n^{10})$$\Rightarrow${$\mathcal{O}(n^2)$} | $\mathcal{O}(n^{10})$$\Rightarrow${$\mathcal{O}(n^3)$} |

**Table 2**. Time for experiments in Fig 3.
|     | SOM      | VIS     | SAL   | DAN    | DMN     | FPN    | LIM | Random |
|-----------|----------|---------|-------|--------|---------|--------|-----|--------|
| 2-way     | 1s       | 2s      | 1s    | 2s     | 2s      | 2s     | 1s  | 12s    |
| 3-way     | 12s      | 18s     | 7     | 13s    | 8s      | 8s     | 1s  | 13s    |
| 4-way     | 5m36s    | 5m      | 2m44s | 2m31s  | 2m      | 2m10s  | 1s  | 2m41s  |
| 5-way     | 2h18m24s | 1h30m9s | 1h3m  | 54m40s | 1h12m4s | 27m16s | 2s  | 29m    |

**Table 3**. Number of terms in the test statistics before and after eliminating the partitions with singletons.
| # of var                 | 4                  | 5                   | 6                    | 7                     | 8                      | 9                        |
|----------------------------------|--------------------|---------------------|----------------------|-----------------------|------------------------|--------------------------|
| Streitberg$\Rightarrow$(centred) | 15$\Rightarrow$(4) | 52$\Rightarrow$(11) | 203$\Rightarrow$(41) | 877$\Rightarrow$(162) | 4140$\Rightarrow$(715) | 21147$\Rightarrow$(3425) |
| Lancaster$\Rightarrow$(centred)  | 12$\Rightarrow$(1) | 27$\Rightarrow$(1)  | 58$\Rightarrow$(1)   | 121$\Rightarrow$(1)   | 248$\Rightarrow$(1)    | 503$\Rightarrow$(1)      |

Meta-review: This paper considers the problem of how to take higher-order interactions into account in statistical modeling, and proposes a hierarchy of d-order interaction measures with kernel-based nonparametric tests to assess their statistical significance. It is also remarkable to draw a connection to the lattice theory.

This paper is overall well-written, and all the reviewers are happy with the contribution of the paper and recommend acceptance. I also agree with it.

I have one comment about related work. This paper has a close relationship to information geometry of hierarchies of probability distributions, which is not discussed in the paper. Information geometry of statistical models with higher-order interactions has been established in a seminal paper [1],  where it is shown that natural parameters of the log-linear model with higher-order interactions and its expectations form mutually orthogonal coordinate systems of (dually flat) statistical manifolds.
Moreover, recently, this line of research has been further generalized into any posets using Zeta/Möbius transformation in Section 4 in [2], where the lattice theory has been already used in probabilistic analysis.
Since they are quite relevant to the contribution of this paper, I strongly recommend adding some discussion in the camera-ready version of the paper. 

To summarize, I would appreciate the contribution of the paper, which is sound and novel. Therefore, I recommend acceptance of the paper.

[1] S. Amari. Information geometry on hierarchy of probability distributions. IEEE Transactions on
Information Theory, 47(5):1701–1711, 2001.  
[2] M. Sugiyama,  et al.: Tensor Balancing on Statistical Manifold, ICML2017.
Decision: Accept (poster)