Exploring the space of graphs with fixed discrete curvatures

Abstract

Discrete curvatures are quantities associated to the nodes and edges of a graph that reflect the local geometry around them. These curvatures have a rich mathematical theory and they have recently found success as a tool to analyze networks across a wide range of domains. In this work, we consider the problem of constructing graphs with a prescribed set of discrete edge curvatures, and explore the space of such graphs. We address this problem in two ways: first, we develop an evolutionary algorithm to sample graphs with discrete curvatures close to a given set. We use this algorithm to explore how other network statistics vary when constrained by the discrete curvatures in the network. Second, we solve the exact reconstruction problem for the specific case of Forman–Ricci curvature. By leveraging the theory of Markov bases, we obtain a finite set of rewiring moves that connects the space of all graphs with a fixed discrete curvature.
 
Keywords: network geometry, discrete curvature, Markov bases, Markov chain Monte Carlo methods, evolutionary algorithm, edge-based network measures

1 Introduction

Networks are used ubiquitously as a tool to model, analyze and design complex systems in a wide range of domains of science and technology such as neuroscience [24], epidemiology [48], economics [28, 6], urban systems [4], energy infrastructure [47] and many more. Across these different disciplines, networks, graphs and their higher-order generalizations [35, 5, 8] provide a unifying mathematical language that allows to translate concepts and tools from one setting to another. Moreover, a series of observations at the turn of the 21st century found that the networks extracted from real-world systems often exhibit common features: heavy-tailed degree distributions [3], a “small world” structure [60] and rich mesoscale connectivity patterns such as communities [26]. This initiated an ongoing pursuit to find models and theories that explain the presence of these features, and develop new mathematical methods to detect, describe and quantify them. We point to [11, 58, 12] for a critical perspective on these observations.

Among the many approaches to study complex networks, this article is concerned with geometry and more specifically with discrete curvature as a tool to describe the structure of graphs; see [10, 62] for an overview of different notions of geometry associated with complex systems and their applications. We now explain in more detail the main topic of the article: discrete curvature.
 
Curvature is a concept from Riemannian geometry, the study of smooth objects such as surfaces equipped with a metric, and quantifies how much these objects differ locally from being flat. While classical definitions of curvature need the structure of a Riemannian manifold, more recent works have generalized these definitions so they can be applied to other metric spaces such as graphs; this results in discrete notions of curvature — see Section 2 for a concrete example.

Definitions of discrete curvature on graphs were given for instance by Forman based on Laplacian operators [23, 53, 32], by Ollivier and later Lin, Lu and Yau based on neighborhood overlap [44, 45, 36], by Bakry and Émery based on the Bochner identity [2], by Erbar and Maas based on a certain entropy functional [21], by Steinerberger based on solutions to an equilibrium problem [55], and by Devriendt and Lambiotte based on the effective resistance [14, 15]. Importantly, many of these discrete curvatures satisfy properties and results analogous to the powerful theorems from differential geometry; see for instance [37, 36, 7, 31, 39, 23, 38, 15, 55]. Of particular relevance to the central question of emergent structure and organization in complex systems, we note that many notions of discrete curvature satisfy powerful ‘local-to-global’ theorems, where local constraints on the curvature lead to statements on the global structure of the graph.

Aside from a rich mathematical theory, discrete curvature has also found successful applications as a tool to analyze and characterize the structure of complex real-world networks. It has been used for instance in the context of community detection [43, 27, 22], in the study of biological networks [51, 19], in characterizing connectivity patterns in brain networks [20], in studying market instabilities in financial networks [49], and to detect and address the problem of oversquashing in graph neural networks [57, 9].
 
In this article, we consider the following question: How to construct graphs with prescribed discrete curvatures? To our knowledge, this inverse problem has not yet been addressed for any notion of discrete curvature, and we consider it a basic methodological gap in trying to understand how the geometry of a graph or network influences its other features. As a first step in addressing this methodological gap, this article provides two different graph reconstruction approaches for specific choices of discrete curvatures. We furthermore explore, both experimentally and theoretically, the resulting space of graphs with fixed discrete curvatures.

In the first part, Section 3, we consider the problem of approximate reconstruction and develop a methodology that works for any notion of discrete curvature. We use a Markov chain Monte Carlo-type algorithm (MCMC) to iteratively construct graphs with discrete curvatures close to a target set of curvatures. We then use this algorithm to explore the ensemble of graphs with a given discrete curvature sequence in a number of experiments. For computational reasons, our discussion focuses on two simple notions of discrete curvature, Forman–Ricci and augmented Forman–Ricci curvature, but we have verified the methodology for other curvatures as well.

In the second part, Section 4, we take a more theoretical approach and consider the problem of exact reconstruction of graphs from a given target Forman–Ricci curvature sequence. In this specific case, the problem reduces to a combinatorial question involving the so-called joint degrees in the graph; see for instance [54]. Using known characterizations of possible joint degrees in combination with the theory of Markov bases [17], we show that there exists a finite set of graph rewiring moves that connect any two graphs with the same Forman–Ricci curvatures and node degrees; this is the content of Theorem 4.12. The set of rewiring moves can be computed using existing computer algebra packages, and we give explicit examples for the case of graphs with a small maximum degree. Our result naturally leads to an ergodic Markov chain over the set of all graphs with equal Forman–Ricci curvatures and degrees.
 
Organization: The rest of this article is organized as follows. Section 2 introduces the relevant definitions from graph theory and discrete curvature and formalizes the reconstruction question. In Section 3, we propose an evolutionary MCMC algorithm for approximate reconstruction from a target curvature sequence. We then use the algorithm in a number of experiments and describe network statistics of graphs sampled using the algorithm. In Section 4, we consider the exact reconstruction problem. We formalize the problem and solve it using the theory of joint degree matrices and Markov bases. Section 5 concludes the article with a summary of the results and an outlook towards future applications and extensions.

2 Background and Forman–Ricci curvatures

We consider simple graphs $G=(V,E)$ where $V$ is the set of nodes and $E$ is the set of edges, which are pairs $\{u,v\}\in E$ of distinct nodes that are connected in the graph. As mentioned in the introduction, we will characterize the structural properties of graphs by considering their discrete curvature. Our focus will be on two definitions of discrete Ricci curvature on edges.

The Forman–Ricci curvature for graphs is defined with respect to a choice of positive weights $w$ for the nodes and edges. The curvature of an edge $e=\{u,v\}$ is equal to [53]

$$F_{w}(e)=w_{u}\Biggl{(}1-\sum_{\begin{subarray}{c}e^{\prime}\ni u\\ e^{\prime}\neq e\end{subarray}}\sqrt{\frac{w_{e}}{w_{e^{\prime}}}}\Biggr{)}+w_{v}\Biggl{(}1-\sum_{\begin{subarray}{c}e^{\prime}\ni v\\ e^{\prime}\neq e\end{subarray}}\sqrt{\frac{w_{e}}{w_{e^{\prime}}}}\Biggr{)},$$

where $e^{\prime}$ runs over all edges that contain $u$ (resp. $v$) and are distinct from $e$. Setting all weights to $w=1$ results in the combinatorial formula for the Forman–Ricci curvature of an edge that will be used in this article:

$$F(e)=4-\deg(u)-\deg(v),$$ (1)

for $e=\{u,v\}$ and where $\deg(u)$ denotes the degree of node $u$, i. e., its number of neighbors in $G$. Clearly, Forman–Ricci curvature is a very simple notion of discrete curvature that only depends on the local degree information of the graph. For this reason, we will also consider an extended notion of discrete Ricci curvature based on Forman’s original curvature for CW complexes [23]. The augmented Forman–Ricci curvature for unweighted graphs is defined in [50] as

$$F^{\#}(e)=4-\deg(u)-\deg(v)+3\mathcal{T}(e),$$ (2)

where $\mathcal{T}(e)$ is the number of triangles in the graph that contain the edge $e=\{u,v\}$. The augmented Forman–Ricci curvature is widely used because it can be computed efficiently while still capturing the local geometry of a graph well; for instance, it correlates highly with other notions of discrete curvature such as Ollivier–Ricci curvature [50] and Bakry–Émery curvature [40]. One can further augment $F^{\#}(e)$ by considering cycles of length four and more, but this involves a tradeoff between expressivity and computational complexity of the discrete curvature; see for instance [32, 56, 22] for more on this. Figure 1 shows a small graph with the Forman–Ricci and augmented Forman–Ricci curvatures computed for all edges.

Figure 1: The Forman–Ricci curvatures $F(e)$ and augmented Forman–Ricci curvatures $F^{\#}(e)$ in a graph. The differences due to the presence of a triangle are highlighted in red.

For a given graph $G$, the curvature sequence is the multiset (set with possible repeats) of all discrete curvatures of its edges. For instance, the Forman–Ricci curvature sequence is the multiset $\{F(e)\}_{e\in E}$, which for the example in Figure 1 equals $\{-3,-3,-2,-2,-1,-1,0,0\}$. We can now make the main methodological question of this article more precise:

Given a multiset $$\mathcal{C}^{\star} \mathcal{C}^{\star}$$ of numbers, construct a graph with curvature sequence . (Q1)

The answer to this question can then be used to address a second goal of the article: to explore the space of all graphs that have a given discrete curvature.

Section 3 addresses question (Q1) algorithmically by proposing an algorithm that samples graphs $\hat{G}$ such that $\mathcal{C}(\hat{G})\approx\mathcal{C}^{\star}$, for any notion of discrete curvature. Section 4 addresses question (Q1) theoretically for the case of Forman–Ricci curvature. We show that there exists a finite set of graph operations, such that any two graphs $G,G^{\prime}$ with the same Forman–Ricci curvature sequence $\mathcal{C}(G)=\mathcal{C}(G^{\prime})$ and degree sequence can be transformed into each other by using these graph operations. Throughout the article, we limit our scope to choices of $\mathcal{C}^{\star}$ for which there exists at least one graph $G$ with curvature sequence $\mathcal{C}^{\star}$. Characterizing the multisets that are realizable as discrete curvature sequences of a graph is an interesting problem in itself, and it is still unexplored for most notions of discrete curvature.

3 Approximate graph reconstruction from curvatures

The main methodological question of this article (Q1) asks to construct graphs with predefined discrete curvatures. Similar to reconstructing graphs from statistics such as motif counts or eigenvalues [29], this is a very difficult task in general and it highly depends on the choice of discrete curvature. Therefore, this section proposes an algorithm for approximate reconstruction that works for any choice of curvature. For a given multiset $\mathcal{C}^{\star}$, the algorithm samples graphs $\hat{G}$ such that $\mathcal{C}(\hat{G})\approx\mathcal{C}^{\star}$. We use this algorithm to explore the ensemble of graphs with curvatures close to a given target sequence. In particular, we study how a number of standard network summary statistics vary over this ensemble.

3.1 Evolutionary algorithm for sampling graphs with given curvatures

We follow the methodology of [29] and use an evolutionary Markov chain Monte Carlo-type algorithm (MCMC) that consists of three steps: initialization, mutation and selection. The algorithm starts with an initial graph $G_{0}$ which is then evolved iteratively to make its curvature sequence gradually closer to the target sequence. After choosing a notion of discrete curvature of interest, the algorithm takes three inputs: (i) the number of nodes $n$ and the target sequence $\mathcal{C}^{\star}$, which is a multiset of numbers, (ii) a parameter $\theta$ which introduces variability to avoid the algorithm from getting stuck and (iii) a parameter $T$, which is the number of steps after which the algorithm halts. In practice we choose the parameters $\theta$ and $T$ for a given sequence using grid search. We now describe the algorithm steps in more detail.
$\blacktriangleright$ Initialization: Construct a random graph $G_{0}$ with $n$ nodes and $m:=|\mathcal{C}^{\star}|$ edges, starting from a random tree (to guarantee connectedness) and with the remaining edges randomly added while avoiding multi-edges. The iteration parameter $t$ is set to zero.
$\blacktriangleright$ Mutation: Select a random edge $\{u,v\}\in E(G_{t})$ in the graph and construct a graph $G^{\prime}$ in which the edge $\{u,v\}$ is rewired to $\{u,v^{\prime}\}$, where $v^{\prime}$ is chosen uniformly at random while avoiding isolated nodes and multi-edges.
$\blacktriangleright$ Selection: Compute the discrete curvature sequences $\mathcal{C}(G_{t})$ and $\mathcal{C}(G^{\prime})$ and compare with the target curvature sequence by computing the mean squared error:

$$\text{MSE}_{t}=\frac{1}{m}\sum_{i=1}^{m}(\mathcal{C}(G_{t})_{i}-\mathcal{C}^{\star}_{i})^{2}\quad\text{~{}and~{}}\quad\text{MSE}^{\prime}=\frac{1}{m}\sum_{i=1}^{m}(\mathcal{C}(G^{\prime})_{i}-\mathcal{C}^{\star}_{i})^{2},$$

where the curvature sequences are sorted in ascending order. The mutation $G^{\prime}$ is accepted if it results in a decrease in the mean squared error $\Delta\text{MSE}=\text{MSE}^{\prime}-\text{MSE}_{t}<0$, or otherwise it is accepted with probability $p=\exp(-\Delta\text{MSE}/(\text{MSE}_{t}\cdot\theta))$. Finally, set $t\leftarrow t+1$ and $G_{t}\leftarrow G^{\prime}$ if the mutation is accepted or $G_{t}\leftarrow G_{t-1}$ otherwise, and continue the process until $t=T$.

Figure 2 illustrates the evolution of the curvature sequences throughout the algorithm. It shows the histogram of the target sequence $\mathcal{C}^{\star}$ and some curvature sequences $\mathcal{C}(G_{t})$ during the algorithm; note the good convergence of the final curvature sequence $\mathcal{C}(G_{T})$. While this methodology is applicable to any choice of discrete curvature, we focus on Forman–Ricci and augmented Forman–Ricci curvature for computational reasons in our experiments. We found that the algorithm still converges for Ollivier–Ricci curvature and resistance curvature target sequences, but with prohibitively long computation times.

Figure 2: Evolution of curvature sequences throughout the sampling algorithm. Top left: Histogram of the target sequence $\mathcal{C}^{\star}$. This is the Forman–Ricci curvature sequence of the C. elegans metabolic network (see Section 3.2.1). Bottom left: Evolution of the mean squared error MSE${}_{t}$ between the target sequence and the Forman–Ricci curvature sequences. Highlight box: Histogram of Forman–Ricci curvatures $\mathcal{C}(G_{t})$ at four points during the algorithm: (A) the initial graph (B) after $5\text{\times}{10}^{4}$ accepted mutations, (C) after $11\text{\times}{10}^{4}$ accepted mutations and (D) at the end of the algorithm, after $16\text{\times}{10}^{4}$ accepted mutations and $T=$2\text{\times}{10}^{6}$$ total iterations.

3.2 Approximate sampling experiments

We now use the algorithm described in Section 3.1 to explore ensembles of graphs with a given discrete curvature sequence. For both Forman–Ricci and augmented Forman–Ricci curvature, we consider four target curvature sequences and sample $1000$ graphs. We then compute a number of network summary statistics and explore how they vary over these graphs.

3.2.1 Target sequences

To guarantee target sequences which are realizable, we use curvature sequences of existing graphs. We consider graphs obtained from three widely-studied random graph models [42] and one real-world network. The considered graphs are as follows:

• Erdős–Rényi (ER) random graph. The graph has $n=500$ nodes and every pair of nodes is connected independently with probability $0.01$. These graphs result in a unimodal target sequence.

• Stochastic block model (SBM) random graph. The graph has $n=500$ nodes divided in two equal sized groups. Pairs of nodes in group $1$ are connected with probability $p=0.05$, in group $2$ with probability $0.01$ and pairs of nodes from different groups with probability $p=0.003$. These graphs result in a bimodal target sequence.

• Barabási–Albert (BA) random graph. The graph has $n=500$ nodes and is constructed by starting from a star graph with four nodes and then adding new nodes. Each new node is connected to three of the existing nodes with probability proportional to their degrees. These graphs result in a heavy-tailed target sequence.

• Real-world network. We consider the metabolic network of the Caenorhabditis elegans worm (C. elegans) [34, 18, 30, 46], this is a well-studied graph in biology. This graph has $n=453$ nodes and $m=2025$ edges and results in a heterogeneous/complex target sequence.

Figure 3: Histograms of the Forman–Ricci curvature of an Erdős–Rényi, stochastic block model and Barabási–Albert random graph and the C. elegans metabolic network.

3.2.2 Summary statistics

To summarize the structure of the graphs obtained from our algorithm, we compute four global network statistics:

• Number of triangles: the number of triangles in a graph is a simple but important summary statistic. It contains information about the density of the graph and it is closely related to transitivity [59] and discrete curvature [31].

• Average clustering coefficient: the clustering coefficient $C_{v}$ of a node $v$ is the number of triangles that contain $v$ divided by the number of possible triangles $\binom{\deg(v)}{2}$. The average clustering coefficient $\frac{1}{n}\sum C_{v}$ is an important graph property that reflects its geometry [31, 33].

• Average shortest path length: this is the average, over all pairs of nodes $u$ and $v$, of the smallest number of steps between $u$ and $v$. It is an intuitive summary of the metric structure of a graph.

• First positive eigenvalue: the normalized graph Laplacian is the $n\times n$ matrix $L=I-D^{-1}A$ where $I$ is the identity matrix, $D$ is the diagonal degree matrix and $A$ is the adjacency matrix. The eigenvalues of the Laplacian contain a lot of information about the graph structure and dynamical processes such as diffusion and random walks associated with the graph [13, 41]. In particular, $L$ is a positive semidefinite matrix and its first positive eigenvalue $\lambda_{2}$ quantifies how well-connected the graph is.

3.2.3 Experiments

Before sampling the ensemble of graphs, we determine the optimal parameters $\theta$ and $T$ for each of the target sequences. Table 1 shows the parameter values obtained via grid search.

		ER	SBM	BA	C. elegans
Forman–Ricci curvature sequence	$\theta$	-.$3\text{\times}{10}^{-3}$	-.$3\text{\times}{10}^{-3}$	-.$3\text{\times}{10}^{-3}$	-.$8\text{\times}{10}^{-4}$
	$T$	$1\text{\times}{10}^{6}$	$1\text{\times}{10}^{6}$	$1\text{\times}{10}^{6}$	$2\text{\times}{10}^{6}$
augmented Forman–Ricci curvature sequence	$\theta$	-.$3\text{\times}{10}^{-3}$	-.$3\text{\times}{10}^{-3}$	-.$4\text{\times}{10}^{-3}$	-.$6\text{\times}{10}^{-4}$
	$T$	$2\text{\times}{10}^{6}$	$2\text{\times}{10}^{6}$	$2\text{\times}{10}^{6}$	$2.5\text{\times}{10}^{6}$

Table 1: Optimal $\theta$ and $T$ parameters determined via grid search.

We now consider the four graphs (ER, SBM, BA, C. elegans), compute their Forman–Ricci curvature sequence as a target sequence $\mathcal{C}^{\star}$ and then compute $1000$ graphs with $\mathcal{C}(\hat{G})\approx\mathcal{C}^{\star}$ using our algorithm. As an indication of the quality of the approximation, Table 2 shows the mean squared errors for the sampled graphs. Generally, we observe that the MSE increases as the graphs become less uniform, as is the case for the Barabási–Albert and the C. elegans network.
 
Figure 4 shows the histograms of the network statistics (in blue) and the value of this statistic for the starting graph (red line). For each statistic, we also plot the standard score $\tfrac{|s^{\star}-\langle s\rangle|}{\sigma}$, where $s^{\star}$ is the value of the statistic for the target graph, $\langle s\rangle$ is the ensemble average of the statistic and $\sigma$ is the standard deviation of the statistic over the ensemble. This number is a proxy for how typical or atypical the statistics of the target graph are with respect to the ensemble.

Figure 4: Network statistics for graphs sampled using the algorithm from Section 3.1, based on four target Forman–Ricci curvature sequences (ER, SBM, BA, C. elegans). The histograms of the observed statistics in the sampled graph ensemble are shown in blue and the red line indicates the statistic of the target graph together with its standard score.

We note a few observations: the standard scores for the ER-based target sequence are all below $1$. This suggests that the ensemble of graphs obtained from the algorithm is not very different from an ensemble of ER random graphs. Second, for the SBM-based target sequence, we see that the average shortest path length is much smaller for the sampled graphs than for the target graph, and that the first positive eigenvalue is much larger. Both observations suggest that the block structure of the SBM graph, i. e., two groups of nodes which are poorly connected, is not well-reproduced by the sampled graphs. Finally, for the C. elegans-based target sequence, we see that all standard scores are very high. This implies that the structure of this graph is very different from the sampled graphs. This is not surprising since only fixing the degree sums of edges is a very local and coarse structural constraint.
 
We repeat the same experiment with the augmented Forman–Ricci curvature. The mean squared error between the sampled and target curvature sequences are shown in Table 2 and the MSE is again higher for the Barabási–Albert graph and the C. elegans network. The results of the experiment are shown in Figure 5.

Figure 5: Network statistics for graphs sampled using the algorithm from Section 3.1, based on four target augmented Forman–Ricci curvature sequences (ER, SBM, BA, C. elegans). The histograms of the observed statistics in the sampled graph ensemble are shown in blue and the red line indicates the statistic of the target graph together with its standard score.

We note a few observations: compared to Figure 4, the standard scores for the ER-based target sequence have all increased. This is not unexpected, since the graphs are sampled based on the augmented Forman–Ricci curvature which takes into account the presence of triangles in the graph, whereas ER graphs have no correlations between different edges. Secondly, the same observation as before holds for the SBM-based target sequence: the smaller shortest path lengths and larger first positive eigenvalue suggest that the two-block structure of the SBM is not reproduced in the sampled graphs.

MSE	Forman				Augmented Forman
	ER	SBM	BA	C. elegans	ER	SBM	BA	C. elegans
min	0.0	0.016	0.038	0.111	0.018	0.153	0.195	15.539
max	0.005	0.039	0.193	0.416	0.191	0.285	3.757	21.469
mean	0.001	0.027	0.088	0.033	0.094	0.212	1.115	18.188

Table 2: Minimal, maximal and mean of the MSE between $\mathcal{C}^{\star}$ and $\mathcal{C}(G_{T})$ for the $1000$ sampled graphs in the experiments of Figures 4 and 5.

Finally, we note that one should take care in interpreting the histograms of network statistics, since the specific design of our algorithm does not give much control or insight into the distribution from which we sample graphs. Furthermore, we do not have any mixing time results on which we can base our choice of stopping parameter $T$. This makes it hard to control the influence of the initialization $G_{0}$ on the ensemble of graphs which we sample from. In future work, these issues could be addressed by designing a more tailored algorithm or by trying to find exact expressions for the probability distribution over graphs, for instance using a soft maximum entropy model [52].

4 Exact graph reconstruction from curvatures and degrees

In this section, we take a second approach to the main question (Q1) and consider the exact reconstruction problem of finding graphs whose discrete curvatures precisely match a given multiset of numbers. This problem clearly depends on the notion of discrete curvature, and we will focus on Forman–Ricci curvature $F(\{u,v\})=4-\deg(u)-\deg(v)$. For reasons which will be explained later, we not only consider fixed curvatures but also fixed degrees. In the main result, Theorem 4.12, we find that there exists a finite set of graph rewiring moves which can transform any two multigraphs with the same Forman–Ricci curvatures and degree sequences into each other. This result follows from encoding the degree and curvature constraints algebraically using a joint degree matrix and using the theory of Markov bases. Our approach then decomposes into two parts: first, finding all joint degree matrices that are compatible with the curvature and degree constraints and second, finding all graphs with the same joint degree matrix.

The rest of this section will use the following conventions: a multigraph can have multi-edges and self-loops while a simple graph cannot, and we write “graph” if the distinction is not relevant. The maximum degree is denoted by $\Delta$, the degree sequence is the $\Delta$-tuple $(\text{\# nodes of degree $a$})_{a=1}^{\Delta}$ and the curvature sequence is the $(2\Delta-1)$-tuple (# edges with Forman–Ricci curvature $4-\kappa$)${}_{\kappa=2}^{2\Delta}$; see Figure 6 for an example. Note that this is different from how we defined the curvature sequence $\mathcal{C}$ before; it records the frequency of all possible curvatures rather than the multiset of curvatures itself. Since both definitions contain the same information, we will use the same name.

4.1 Joint degree matrix and transpositions

The Forman–Ricci curvature of an edge depends on the sum of the degrees of its end points. This information can be represented by the following matrix: the joint degree matrix (JDM) of a graph $G$ is a symmetric $\Delta\times\Delta$ matrix $J$ with entries

$$J_{ab}=J_{ba} =\left|\big{\{}\{u,v\}\in E(G)\,:\,\{\deg(u),\deg(v)\}=\{a,b\}\big{\}}\right|$$

$$=\text{~{}\# edges whose end points have degrees $a$ and $b$.}$$

We note that the degree and curvature sequences of a graph can be computed from its JDM: the number of nodes with degree $a$ is equal to the row sum $(\sum_{b}J_{ab}+J_{aa})/a$ and the number of edges with Forman–Ricci curvature $\kappa$ is equal to the anti-diagonal sum $\sum_{a\leq b,a+b=\kappa}J_{ab}$. Figure 6 shows the JDM and associated sequences for a small graph that will be our running example.

Figure 6: A graph $G$ and its joint degree matrix (JDM). The figures on the right show the degree and curvature sequence of $G$ and highlight how these can be obtained from the JDM.

The following result from [54] gives a complete characterization of JDMs.

Theorem 4.1.

A symmetric matrix $J\in\mathbb{N}^{\Delta\times\Delta}$ is the JDM of a simple graph if and only if it satisfies, for all distinct $a,b\in\{1,\dots,\Delta\}$, the conditions

• $|V_{a}|=(\sum_{c=1}^{\Delta}J_{ac}+J_{aa})/a$   is an integer,

• $J_{ab}\leq|V_{a}|\cdot|V_{b}|$, and

• $J_{aa}\leq\binom{|V_{a}|}{2}$.

The authors of [54] gave an algorithm to construct a simple graph from a given JDM. More importantly, they also showed that the set of all graphs with the same JDM is connected through the graph operation shown in Figure 7: let $\mathcal{E}=\big{\{}\{u,v\},\{x,y\}\big{\}}\subseteq E(G)$ be two edges such that $\deg(u)=\deg(x)$ and let $\mathcal{E}^{\prime}=\big{\{}\{u,y\},\{x,v\}\big{\}}$. The transposition of the edges $\mathcal{E}$ in $G$ is a new graph $G^{\prime}$ obtained by rewiring these two edges. Its nodes and edges are

$$V(G^{\prime})=V(G)\quad\text{~{}and~{}}\quad E(G^{\prime})=\big{(}E(G)\backslash\mathcal{E}\big{)}\cup\mathcal{E}^{\prime}.$$

Figure 7: A transposition rewires two edges in a graph without changing the joint degrees. In this example, edges $\{u,v\},\{x,y\}$ are rewired to $\{u,y\},\{x,v\}$.

Transpositions do not change the JDM of a graph and the theory of matchings in convex bipartite graphs guarantees that the set of multigraphs with the same JDM is connected through transpositions [16, 54]. In general, transpositions in simple graphs may create multi-edges or self-loops depending on the choice of $\mathcal{E}$, but the following result shows that there always exists a choice for which this is not the case.

Theorem 4.2 ([54]).

The set of all simple graphs with a given JDM is connected through transpositions.

In other words, all simple graphs with the same JDM can be found by performing transpositions, in such a way that every intermediate step is also a simple graph. Since transpositions do not change the JDM, these graphs will also have the same curvature and degree sequences. This is the first step towards our main result.
 
We note the following count for the number of multigraphs with a given JDM:

Proposition 4.3.

The number of node-labeled multigraphs with a given JDM is equal to

$$\prod_{1\leq a\leq\Delta}\frac{(a|V_{a}|)!}{(a!)^{|V_{a}|}\cdot 2^{J_{aa}}}\cdot\prod_{1\leq a\leq b\leq\Delta}\frac{1}{J_{ab}!}.$$ (3)

Proof. We have $|V_{a}|$ nodes of degree $a$. Each such node has $a$ “stubs” that need to be matched to $a$ of the $a|V_{a}|$ “half-edges”. For the first node, there are $\binom{a|V_{a}|}{a}$ ways to match its stubs, for the second node there are $\binom{a(|V_{a}|-1)}{a}$ ways, and so on. This gives a total of

$$\binom{a|V_{a}|}{a}\binom{a(|V_{a}|-1)}{a}\dots\binom{a}{a}=\frac{(a|V_{a}|)!}{a!(a(|V_{a}|-1))!}\frac{(a(|V_{a}|-1)!}{a!(a(|V_{a}|-2))!}\dots\frac{a!}{a!}=\frac{(a|V_{a}|)!}{(a!)^{|V_{a}|}}$$

ways to connect nodes of degree $a$ to half-edges. This counting procedure distinguishes the half-edges, so to obtain edge-unlabeled multigraphs we have to forget the half-edge labels. This makes any permutation of the $J_{ab}$ edges which are connected to nodes with degree $a$ and $b$ indistinguishable, so we divide the count by $J_{ab}!$. Second, forgetting the half-edge labels makes any flip of the $J_{aa}$ edges connected to two nodes of equal degree $a$ indistinguishable, so we divide the count by $2^{J_{aa}}$. This results in (3) and completes the proof. ∎

For the JDM in Figure 6 with degree sequence $(2,2,2,1)$, the number (3) equals $\frac{2!4!6!4!}{2^{5}(3!)^{2}4!}=30$, so there are $30$ node-labeled multigraphs with this JDM. In practice, we may be interested in simple and non-isomorphic graphs for which this number serves as an upper bound.

4.2 Markov bases

In the previous section, we found that all graphs with a given JDM are connected by transpositions. As a second step, we now want to find all JDMs which are compatible with a given curvature and degree sequence. More precisely, these are all nonnegative integer matrices $J$ that satisfy conditions (i)—(iii) in Theorem 4.1 (valid JDMs) and have prescribed row sums (degrees) and anti-diagonal sums (curvatures). We will show how this set of matrices can be obtained using Markov bases.

The theory of Markov bases was developed by Diaconis and Sturmfels in [17] to solve the problem of sampling contingency tables, i. e., nonnegative integer matrices with fixed row and column sums, in the context of statistics. The main result in this theory is that this sampling problem can be solved using concepts and tools from commutative algebra. We first give the main definition and theorem and then apply it to our problem.

Definition 4.4.

Let $A\in\mathbb{Z}^{d\times n}$ be an integer matrix. A Markov basis for $A$ is a finite set $\mathcal{M}=\{\beta_{1},\dots,\beta_{M}\}\subseteq\ker_{\mathbb{Z}}(A)$ of integer vectors in the kernel of $A$, such that for every $u,v\in\mathbb{N}^{n}$ with $Au=Av$, there exists a sequence $(i_{1},s_{1}),(i_{2},s_{2}),\dots,(i_{L},s_{L})\in\{1,\dots,M\}\times\{\pm 1\}$ that satisfies

• $u+\sum_{j=1}^{L}s_{j}\beta_{i_{j}}=v$, and

• $u+\sum_{j=1}^{\ell}s_{j}\beta_{i_{j}}\geq 0\text{~{}for all~{}}1<\ell<L$.

For $u\in\mathbb{N}^{n}$ we call $\mathcal{F}_{A}(u)=\{v\in\mathbb{N}^{n}:Au=Av\}$ the fiber of $u$. The adjective ‘Markov’ reflects that a Markov basis can be used to design a Markov chain on a fiber: points in the fiber are states of the Markov chain and the transitions out of any state $x$ are those vectors $\beta$ in the Markov basis for which $x+\beta\geq 0$. The Markov basis properties guarantee that this chain is ergodic.

The question of how to find a Markov basis for a given integer matrix is addressed by the fundamental theorem of Markov bases (see also Appendix A).

Theorem 4.5 (Fundamental theorem of Markov bases [17]).

Let $A\in\mathbb{Z}^{d\times n}$, then $\mathcal{M}\subseteq\ker_{\mathbb{Z}}(A)$ is a Markov basis for $A$ if and only if the binomials $\{\mathbf{x}^{\beta^{+}}-\mathbf{x}^{\beta^{-}}:\beta\in\mathcal{M}\}$ generate the ideal $\langle\mathbf{x}^{\beta^{+}}-\mathbf{x}^{\beta^{-}}:\beta\in\ker_{\mathbb{Z}}(A)\rangle$.

Appendix A provides further background on the relevant algebraic concepts, but the details of the theorem are not relevant for the purpose of this article. The most important fact is that a finite Markov basis always exists and that it can be computed using computer algebra systems for sufficiently small matrices $A$. For instance, using the software package 4ti2 [1] we compute that

$$A=\begin{pmatrix}1&1&1&1\\ 0&1&2&3\end{pmatrix}\text{~{}has Markov basis~{}}\mathcal{M}=\left\{\begin{pmatrix}[r]0\\ 1\\ -2\\ 1\end{pmatrix},\begin{pmatrix}[r]1\\ -2\\ 1\\ 0\end{pmatrix},\begin{pmatrix}[r]1\\ -1\\ -1\\ 1\end{pmatrix}\right\}.$$

For larger matrices, these computations can become infeasible in practice and other approaches such as relaxed notions of Markov bases or restrictions to particular fibers may be necessary [25].
 
We now show that the set of all JDMs with a given degree and curvature sequence is the fiber for some matrix $A_{\Delta}$. This means that the associated Markov basis $\mathcal{M}_{\Delta}$ can be used to explore this fiber. We stress that this matrix and its Markov basis only depend on the maximum degree $\Delta$, and not on the specific degree and curvature sequences.

We start by encoding a JDM and its constraints as a vector $\gamma\in\mathbb{N}^{n}$, where $n=2\binom{\Delta+1}{2}$. This vector contains the upper-triangular JDM entries $(J_{ab})_{a\leq b}$ and the slack variables $(S_{ab})_{a\leq b}$ which we define as

$$S_{aa}=\binom{|V_{a}|}{2}-J_{aa}\quad\text{~{}and~{}}\quad S_{ab}=|V_{a}|\cdot|V_{b}|-J_{ab}\quad\text{~{}for $1\leq a<b\leq\Delta$,}$$

with $|V_{a}|$ defined as before from $J$. The slack variables are used to encode the inequalities on the JDM entries given in Theorem 4.1: a matrix $J$ is a valid JDM if and only if the associated slack variables are nonnegative integers, and thus $\gamma\in\mathbb{N}^{n}$. The matrix $J$ can be reconstructed from $\gamma$ by forgetting the slack variables, so we can write $J(\gamma)$.

Next, to construct the matrix $A_{\Delta}$, we consider the equations in the $J$ and $S$ variables that are invariant when fixing the degree and curvature sequences:

• $\sum_{b}J_{ab}+J_{aa}$ for all $a\in\{1,\dots,\Delta\}$ (fixed degree sequence),

• $\sum_{a\leq b,a+b=\kappa}J_{ab}$ for all $\kappa=\{2,\dots,2\Delta\}$ (fixed curvature sequence),

• $S_{ab}+J_{ab}$ for all $1\leq a\leq b\leq\Delta$ (definition of slack variables).

These $d=\binom{\Delta+1}{2}+3\Delta-1$ equations are linear functions of $\gamma$ with integer coefficients and they determine the rows of the $d\times n$ matrix $A_{\Delta}$.

To illustrate, we give $A_{\Delta}$ explicitly for $\Delta=4,5$ in Examples 4.8 and 4.9 below. The relevance of this matrix is shown by the following proposition.

Proposition 4.6.

Let $J$ be a JDM with corresponding vector $\gamma$. The fiber $\mathcal{F}_{A_{\Delta}}(\gamma)$ contains all JDMs with the same degree and curvature sequences as $J$.

Proof. Let $J$ be a JDM with corresponding vector $\gamma$. We consider the fiber

$$\mathcal{F}_{A_{\Delta}}(\gamma)=\{\gamma^{\prime}\in\mathbb{N}^{n}\,:\,A_{\Delta}\gamma^{\prime}=A_{\Delta}\gamma\}.$$

Let $\gamma^{\prime}\in\mathcal{F}_{A_{\Delta}}$ and write $J^{\prime}_{ij}$ and $S^{\prime}_{ij}$ for the corresponding entries in this vector. From the rows of $A_{\Delta}$ that encode the fixed degree sequence, we obtain

$$a|V^{\prime}_{a}|:=\sum_{b}J^{\prime}_{ab}+J^{\prime}_{aa}=\sum_{b}J_{ab}+J_{aa}=a|V_{a}|\text{~{}for all $a\in\{1,\dots,\Delta\}$}.$$ (4)

From the rows of $A_{\Delta}$ corresponding to the off-diagonal slack equations ($1\leq a<b\leq\Delta)$ and the fact that $\gamma^{\prime}\in\mathbb{N}^{n}\geq 0$ and thus $S^{\prime}_{ab}\geq 0$, we know that

$$S^{\prime}_{ab}+J^{\prime}_{ab}=S_{ab}+J_{ab}=|V_{a}|\cdot|V_{b}|=|V^{\prime}_{a}|\cdot|V^{\prime}_{b}|\Longrightarrow J^{\prime}_{ab}\leq|V^{\prime}_{a}|\cdot|V^{\prime}_{b}|.$$

Similarly, from the diagonal slack equations ($1\leq a\leq\Delta$) and $S^{\prime}_{aa}\geq 0$, we obtain

$$S^{\prime}_{aa}+J^{\prime}_{aa}=S_{aa}+J_{aa}=\binom{|V_{a}|}{2}=\binom{|V^{\prime}_{a}|}{2}\Longrightarrow J^{\prime}_{aa}\leq\binom{|V^{\prime}_{a}|}{2}.$$

By Theorem 4.1, this proves that the matrix $J^{\prime}=J(\gamma^{\prime})$ is a JDM. The degree sequences $(|V^{\prime}_{a}|)_{1\leq a\leq\Delta}$ and $(|V_{a}|)_{1\leq a\leq\Delta}$ are equal by (4), and the rows of $A_{\Delta}$ which encode the fixed curvature sequence guarantee that

$$\sum_{\begin{subarray}{c}{a+b}=\kappa\\ a\leq b\end{subarray}}J^{\prime}_{ab}=\sum_{\begin{subarray}{c}{a+b}=\kappa\\ a\leq b\end{subarray}}J_{ab}$$

and thus that the curvature sequences of $J$ and $J^{\prime}$ are equal. This holds for every $\gamma^{\prime}\in\mathcal{F}_{A_{\Delta}}(\gamma)$ which proves that every vector in the fiber indeed corresponds to a JDM with degree and curvature sequences equal to $J$. The same approach shows that every such JDM is also an element of the fiber, which completes the proof. ∎

Following Proposition 4.6, we can thus use a Markov basis $\mathcal{M}_{\Delta}$ corresponding to the matrix $A_{\Delta}$ to obtain all JDMs with a fixed curvature and degree sequence. At this point, it has become clear why we need to consider fixed degrees: to encode the inequalities in the characterization of JDM matrices using slack variables as $S_{ab}+J_{ab}$, we need $\binom{|V_{a}|}{2}$ and $|V_{a}|\cdot|V_{b}|$ to be constant.

A first result on the Markov basis $\mathcal{M}_{\Delta}$ follows from the specific structure of $A_{\Delta}$.

Proposition 4.7.

The entries of a vector $\beta\in\mathcal{M}_{\Delta}$ corresponding to the slack variables are the negative of the entries corresponding to the $J$ variables.

Proof. By construction, the equation $J_{ab}+S_{ab}$ is invariant so any change $J_{ab}+x$ when applying $\beta$ must be compensated by the opposite change $S_{ab}-x$. ∎

Following Proposition 4.7, we can just consider the $J$-entries of vectors in the Markov basis and we will furthermore rearrange these as a symmetric $\Delta\times\Delta$ matrix and write $\delta J\in\mathcal{M}_{\Delta}$. We now consider the cases $\Delta=4,5$ in more detail to illustrate the definitions and results.

Example 4.8 ($\Delta=4$).

The matrix $A_{4}$ is a $21\times 20$ matrix. It can be written as the block-matrix

$$A_{4}=\begin{pmatrix}B&0\\ I_{10}&I_{10}\end{pmatrix},$$

where $I_{10}$ is the $10\times 10$ identity matrix, $0$ is the zero matrix and $B$ is the $11\times 10$ matrix whose rows encode the fixed degree sequences (first $4$ rows) and fixed curvature sequences (last $7$ rows) in terms of the $J$ variables:

$$\setstacktabbedgap{.4em}\setstackgap{L}{.9em}\TABstackTextstyle{\scriptsize\color[rgb]{1,0,0}}\fixTABwidth{T}\savestack\collabels{\tabbedCenterstack{11&121314222324333444}}\savestack\rowlabels{\tabbedCenterstack[c]{|V_{1}|\\ |V_{2}|\\ |V_{3}|\\ |V_{4}|\\ \kappa=2\\ \kappa=3\\ \kappa=4\\ \kappa=5\\ \kappa=6\\ \kappa=7\\ \kappa=8}}\ensurestackMath{B=\stackon{\parenMatrixstack{\makebox[0.0pt]{$2$}&111\\ 1211\\ 1121\\ 1112\\ 1\\ 1\\ 11\\ 11\\ 11\\ 1\\ 1}}{\collabels}\rowlabels}$$

The top indices show the correspondence between the columns of $A_{4}$ and the entries of the $J$ matrix, while the indices on the right indicate the equation encoded by the row. For instance, the second row reflects the sum $2J_{22}+J_{12}+J_{23}+J_{24}$; this equation fixes $2\cdot|V_{1}|$ and thus the number of nodes with degree $a=2$. The $7$th row reflects the sum $J_{13}+J_{22}$; this equation fixes the number of edges with Forman–Ricci curvature $4-\kappa=4-4=0$.

The matrix $A_{4}$ has rank $19$ (see also Proposition 4.10) and a one-dimensional kernel. The corresponding Markov basis $\mathcal{M}_{4}$ thus consists of a single vector

$$\mathcal{M}_{4}=(0,0,1,-1,-1,1,1,-1,0,0)^{T},$$

where we only show the entries corresponding to the $J$ variables, following Proposition 4.7. We can also represent this vector as the symmetric matrix

$$\delta J=\begin{pmatrix}[r]0&0&1&-1\\ 0&-1&1&1\\ 1&1&-1&0\\ -1&1&0&0\end{pmatrix}.$$

Example 4.9 ($\Delta=5$).

The matrix $A_{5}$ is a $29\times 30$ matrix. It can be written in block-matrix form as

$$A_{5}=\begin{pmatrix}B&0\\ I_{15}&I_{15}\end{pmatrix},$$

where $I_{15}$ is the $15\times 15$ identity matrix, $0$ is the zero matrix of appropriate dimensions and $B$ is the $14\times 15$ matrix whose rows encode the fixed degree sequences (first $5$ rows) and fixed curvature sequences (last $9$ rows) in terms of the $J$ variables:

$$\setstacktabbedgap{.4em}\setstackgap{L}{.9em}\TABstackTextstyle{\scriptsize\color[rgb]{1,0,0}}\fixTABwidth{T}\savestack\collabels{\tabbedCenterstack{11&1213141522232425333435444555}}\savestack\rowlabels{\tabbedCenterstack[c]{|V_{1}|\\ |V_{2}|\\ |V_{3}|\\ |V_{4}|\\ |V_{5}|\\ \kappa=2\\ \kappa=3\\ \kappa=4\\ \kappa=5\\ \kappa=6\\ \kappa=7\\ \kappa=8\\ \kappa=9\\ \kappa=10}}\ensurestackMath{B=\stackon{\parenMatrixstack{\makebox[0.0pt]{$2$}&1111\\ 12111\\ 11211\\ 1121\\ 11112\\ 1\\ 1\\ 11\\ 11\\ 111\\ 11\\ 11\\ 1\\ 1}}{\collabels}\rowlabels}$$

For instance, the first row encodes the equation $2J_{11}+J_{12}+J_{13}+J_{14}+J_{15}$, which keeps the number $|V_{1}|$ of nodes with degree $a=1$ fixed. The $8$th row encodes the equation $J_{13}+J_{22}$, which keeps the number of edges with Forman–Ricci curvature $4-\kappa=4-4=0$ fixed.

The matrix $A_{5}$ has rank 27 and using 4ti2 we compute the following Markov basis $\mathcal{M}_{5}$ consisting of $9$ elements:

$$\delta J_{1} = \begin{pmatrix}[r]0&0&0&0&0\\ \phantom{-}0&0&0&1&-1\\ 0&0&-1&1&1\\ 0&1&1&-1&0\\ 0&-1&1&0&0\end{pmatrix},\,\qquad \delta J_{2} = \begin{pmatrix}[r]0&0&0&1&-1\\ 0&0&-1&0&1\\ 0&-1&1&-1&0\\ 1&0&-1&0&0\\ -1&1&0&0&0\end{pmatrix},\,\qquad \delta J_{3} = \begin{pmatrix}[r]0&0&0&1&-1\\ 0&0&-1&1&0\\ 0&-1&0&0&1\\ 1&1&0&-1&0\\ -1&0&1&0&0\end{pmatrix},\,\qquad$$

$$\delta J_{4} = \begin{pmatrix}[r]0&0&1&-2&1\\ 0&-1&2&0&0\\ 1&2&-1&0&-1\\ -2&0&0&1&0\\ 1&0&-1&0&0\end{pmatrix},\,\qquad \delta J_{5} = \begin{pmatrix}[r]0&0&1&-1&0\\ 0&-1&1&0&1\\ 1&1&0&-1&-1\\ -1&0&-1&1&0\\ 0&1&-1&0&0\end{pmatrix},\,\qquad \delta J_{6} = \begin{pmatrix}[r]0&0&1&-1&\phantom{-}0\\ 0&-1&1&1&0\\ 1&1&-1&0&0\\ -1&1&0&0&0\\ 0&0&0&0&0\end{pmatrix},\,\qquad$$

$$\delta J_{7} = \begin{pmatrix}[r]0&0&1&0&-1\\ 0&-1&0&0&2\\ 1&0&1&-2&-1\\ 0&0&-2&1&0\\ -1&2&-1&0&0\end{pmatrix},\,\qquad \delta J_{8} = \begin{pmatrix}[r]0&0&1&0&-1\\ 0&-1&0&1&1\\ 1&0&0&-1&0\\ 0&1&-1&0&0\\ -1&1&0&0&0\end{pmatrix},\,\qquad \delta J_{9} = \begin{pmatrix}[r]0&0&1&0&-1\\ 0&-1&0&2&0\\ 1&0&-1&0&1\\ 0&2&0&-1&0\\ -1&0&1&0&0\end{pmatrix}.\,\qquad$$

If we define $\deg(\delta J):=\frac{1}{2}\sum_{a\leq b}|J_{ab}|$, we see that there are five bases with $\deg(\delta J)=3$, two bases with $\deg(\delta J)=4$ and two bases with $\deg(\delta J)=5$. We will see later that this corresponds to the number of edges involved in the graph operation associated with $\delta J$.

We now give some further results related to the matrix $A_{\Delta}$.

Proposition 4.10.

$\operatorname{rank}(A_{\Delta})=\binom{\Delta+1}{2}+3\Delta-3=d-2$ for all $\Delta\geq 3$.

Proof. We first give two independent linear dependencies between the rows of $A_{\Delta}$ to show that $\operatorname{rank}(A_{\Delta})\leq d-2$. The first linear relation between the rows is

$$\sum_{a=1}^{\Delta}(\text{row $|V_{a}|$})-2\sum_{\kappa=2}^{2\Delta}(\text{row $\kappa$})=\sum_{a=1}^{\Delta}\Biggl{(}\sum_{c}J_{ac}+J_{aa}\Biggr{)}-2\sum_{\kappa=2}^{2\Delta}\Biggl{(}\sum_{\begin{subarray}{c}a+b=\kappa\\ a\leq b\end{subarray}}J_{ab}\Biggr{)}=0.$$

The second linear relation between the rows is

$$\sum_{a=1}^{\Delta}a\times(\text{row $|V_{a}|$})-\sum_{\kappa=2}^{2\Delta}\kappa\times(\text{row $\kappa$})=\sum_{a=1}^{\Delta}a\left(\sum_{c}J_{ac}+J_{aa}\right)-\sum_{\kappa=2}^{2\Delta}\kappa\Biggl{(}\sum_{\begin{subarray}{c}a+b=\kappa\\ a\leq b\end{subarray}}J_{ab}\Biggr{)}=0.$$

To prove that the latter relation holds, we simplify the first sum to $\sum_{v\in V}\deg(v)^{2}$ and the second sum to $\sum_{e\in E}\kappa(e)$, where $\kappa(\{u,v\})=\deg(u)+\deg(v)$. We further simplify

$$\sum_{e\in E}\kappa(e)=\sum_{\{u,v\}\in E}(\deg(u)+\deg(v))=\sum_{v\in V}\sum_{u\sim v}\deg(v)=\sum_{v\in V}\deg(v)^{2},$$

which confirms that the second linear relation equals zero. Next, we show that the matrix $A_{\Delta}$ has a nonsingular $(d-2)$-minor $\tilde{A}$. This implies that $\operatorname{rank}(A_{\Delta})\geq d-2$ which then completes the proof.

First, note that $A_{\Delta}$ can be written as $A_{\Delta}=\left(\begin{smallmatrix}B&0\\ I&I\end{smallmatrix}\right)$, where $0$ is the zero matrix, $I$ is the identity matrix of size $\binom{\Delta+1}{2}$ and $B$ a $(\Delta-1)\times\binom{\Delta+1}{2}$ matrix. The rows corresponding to the identity matrices encode the fixed slack variables and the rows corresponding to $B$ encode the fixed degree and curvature sequences; as a result, the columns of $B$ correspond to $J$ variables. See Examples 4.8 and 4.9 for a concrete example. We show that $B$ has a nonsingular $(3\Delta-3)$-minor $\tilde{B}$ which implies that $A$ has the required nonsingular $(d-2)$-minor $\tilde{A}=\left(\begin{smallmatrix}\tilde{B}&0\\ I&I\end{smallmatrix}\right)$.

Consider the following submatrices of $B$:

$$\setstacktabbedgap{-0.2em}\setstackgap{L}{.9em}\TABstackTextstyle{\scriptsize\color[rgb]{1,0,0}}\fixTABwidth{T}\savestack\collabels{\tabbedCenterstack{\Delta\Delta&(\Delta-1)\Delta(\Delta-1)(\Delta-1)}}\savestack\rowlabels{\tabbedCenterstack[l]{\kappa=2\Delta\\ \kappa=2\Delta-1\\ \kappa=2\Delta-2}}\ensurestackMath{B_{1}=\stackon{\parenMatrixstack{\makebox[0.0pt]{1}&00\\ 010\\ 001}}{\collabels}\rowlabels}$$

$$\setstacktabbedgap{0.4em}\setstackgap{L}{.9em}\TABstackTextstyle{\scriptsize\color[rgb]{1,0,0}}\fixTABwidth{T}\savestack\collabels{\tabbedCenterstack{(\Delta-i)(\Delta+2-i)&(\Delta-i)(\Delta+1-i)(\Delta-i)(\Delta-i)}}\savestack\rowlabels{\tabbedCenterstack[l]{|V_{\Delta-i}|\\ \kappa=2\Delta-2i+1\\ \kappa=2\Delta-2i}}\ensurestackMath{B^{i}=\stackon{\parenMatrixstack{\makebox[0.0pt]{1}&12\\ 010\\ 001}}{\collabels}\rowlabels}$$

for $i\in\{2,\dots,\Delta-1\}$. After reordering rows and columns, we can write matrix $\tilde{B}$ as

$$\tilde{B}=\left(\begin{array}[]{ccccc}B_{1}\\ &B^{2}&&\text{\huge$\star$}\\ &&\ddots\\ &\text{\huge 0}&&B^{(\Delta-1)}\end{array}\right),$$

where $\star$ indicates entries that are not relevant. This is by construction a ($3\Delta-3$)-minor of $B$. It is an upper-triangular matrix with ones along the diagonal, which means that it is nonsingular. As noted, this implies that the $(d-2)$-minor $\tilde{A}$ is nonsingular, which completes the proof. ∎

Understanding the size of the fiber $|\mathcal{F}_{A_{\Delta}}(\gamma)|$ is a hard problem in general. However, making use of the special property that $\ker(A_{4})$ is one-dimensional, we obtain the following bound on the size of the fiber in this case.

Proposition 4.11.

Let $J$ be a JDM with $\Delta=4$. Then there are at most

$$|\mathcal{F}_{A_{4}}(\gamma(J))|\leq 1+\min\{J_{22},J_{14},J_{33}\}+\min\{J_{13},J_{23},J_{24}\}$$

JDMs, including $J$, with the same degree and curvature sequences as $J$.

Proof. We recall from Example 4.8 above that the Markov basis for $A_{4}$ contains a single element

$$\delta J=\begin{pmatrix}[r]0&0&1&-1\\ 0&-1&1&1\\ 1&1&-1&0\\ -1&1&0&0\end{pmatrix}.$$

This means that every JDM with the same degree and curvature sequences as $J$ is of the form $J_{k}=J+k\times\delta J$ with $k\in\mathbb{Z}$. For $J_{k}$ to be a JDM, it is necessary that $J_{k}\geq 0$ and for this reason we know that $k\leq\min\{J_{22},J_{14},J_{33}\}$ and $k\geq\min\{J_{13},J_{23},J_{24}\}$. These conditions are necessary but in general not sufficient, since $J_{k}\geq 0$ alone does not guarantee that the slack variables are also nonnegative. The possible JDMs from the proposition are thus a subset of $J_{k}$ for $k$ within the given bounds, which completes the proof. ∎

We note that the bound is the same for all elements of the fiber. To illustrate that the bound can be strict, consider the JDM in the running example in Figure 6. Here, $J+\delta J$ is a valid JDM but $J-\delta J$ is not and we thus have just two JDMs in the fiber, whereas $1+\min\{J_{22},J_{14},J_{33}\}+\min\{J_{13},J_{23},J_{24}\}=3$.

4.3 Markov moves and main result

To combine Theorem 4.2 on transpositions with Proposition 4.6 on JDMs and Markov bases, we make the observation that applying an element of the Markov basis $J\rightarrow J+\delta J$ to a JDM can be interpreted as a graph operation; this is illustrated in Figure 8.

Figure 8: The figure left illustrates how a Markov basis element $\delta J\in\mathcal{M}_{\Delta}$ can be interpreted as a graph operation that rewires certain edges; we call this a Markov move. The figure right shows a Markov move on the running example.

We define a Markov move as follows: Let $G$ be a graph and $\delta J$ an element of the Markov basis such that $J(G)+\delta J$ is a JDM. Consider a pair $(\mathcal{E},\pi)$ where $\mathcal{E}=(\{i_{\ell},j_{\ell}\})_{\ell=1}^{\deg(\delta J)}\subseteq E(G)$ is a subset of edges that has $|\delta J_{ab}|$ edges with degrees $a,b$ whenever $\delta J_{ab}<0$, and where $\pi$ is a permutation such that the set $\mathcal{E}^{\prime}=(\{i_{\ell},\pi(j_{\ell})\})_{\ell=1}^{\deg(\delta J)}$ has $|\delta J_{ab}|$ edges with degrees $a,b$ whenever $\delta J_{ab}>0$. Existence of such a pair follows from the definition of the Markov basis $\delta J$. A Markov move of $G$ is a new graph $G^{\prime}$ with nodes and edges

$$V(G^{\prime})=V(G)\quad\text{~{}and~{}}\quad E(G^{\prime})=(E(G)\backslash\mathcal{E})\cup\mathcal{E}^{\prime}.$$

See Figure 8 for an example. We now know that all multigraphs with a given degree and curvature sequence are connected by a finite set of graph operations: transpositions, which leave the JDM unchanged, and Markov moves corresponding to the elements of $\mathcal{M}_{\Delta}$, which change the JDM. This is summarized in the following theorem.

Theorem 4.12.

The set of multigraphs with a given Forman–Ricci curvature and degree sequence is connected through transpositions and Markov bases. For small $\Delta$, there exist Markov bases of size

$\Delta$	1	2	3	4	5	6	7	8	$\geq 9$
$|\mathcal{M}_{A_{\Delta}}|$	0	0	0	1	9	111	2662	171964	?

Proof. Let $G$ and $G^{\prime}$ be two multigraphs with the same curvature and degree sequences. If $G$ and $G^{\prime}$ have the same JDM then Theorem 4.12 in the case of multigraphs (see also [54, 16]) says that $G$ can be transformed into $G^{\prime}$ by transpositions. If $G$ and $G^{\prime}$ do not have the same JDM, i. e., $J\neq J^{\prime}$, then by Proposition 4.6 and the fundamental theorem of Markov bases, matrix $J$ can be transformed into $J^{\prime}$ by adding or subtracting a sequence of elements from the Markov basis $\mathcal{M}_{\Delta}$. This implies that $G$ can be transformed using Markov moves into a graph $G^{\prime\prime}$ such that $J^{\prime\prime}=J$. Using Theorem 4.12 again, $G^{\prime\prime}$ can then be transformed into $G^{\prime}$ using transpositions. All intermediate steps in this procedure can be guaranteed to be multigraphs. This concludes the proof of the first part. The results on the size of the Markov bases follows from computations: for $\Delta\in\{1,2,3\}$ the matrix $A_{\Delta}$ has full rank and for $\Delta\in\{4,5,6,7,8\}$ we computed the Markov bases using the software 4ti2. ∎

We note that the set of simple graphs with a given curvature and degree sequence sits inside this bigger set of multigraphs. So in particular, the set of simple graphs is connected by transpositions and Markov moves, with intermediate steps potentially being multigraphs. It is an interesting open question whether all intermediate steps can be simple graphs as well, which would imply a version of Theorem 4.12 for simple graphs.

Finally, we note that Theorem 4.12 may be used to sample graphs with a given degree and curvature sequence. Algorithm 1 shows one possible approach.

Algorithm 1 Sampling graphs with given curvature and degree sequence

graph $G$, Markov basis $\mathcal{M}_{\Delta}$, number of steps $T$, ratio of Markov moves to transpositions $\Theta\in[0,1]$

$t\leftarrow 0$

$J\leftarrow$ JDM of $G$

$S\leftarrow$ slack variables

while $t\leq T$ do

     $p\leftarrow$ sample random variable in $[0,1]$

     if $p\leq\Theta$ then $\triangleright$ do a Markov move

         $\delta J\leftarrow$ random element in $\pm\mathcal{M}_{\Delta}$

         while $J+\delta J\not\geq 0$ or $S-\delta J\not\leq 0$ do

              $\delta J\leftarrow$ random element in $\pm\mathcal{M}_{\Delta}$

         end while

         $G^{\prime}\leftarrow$ Markov move of $G$ associated with $\delta J$

         $(G,J,S)\leftarrow(G^{\prime},J+\delta J,S-\delta J)$

     else if $p>\Theta$ then $\triangleright$ do a transposition

         $G\leftarrow$ transposition of $G$

     end if

     $t\leftarrow t+1$

end while

Return: graph $G^{\prime}$

The performance of any sampling algorithm based on transpositions and Markov moves will depend on how fast the associated Markov chain mixes. The mixing behavior of transpositions is partially understood via its relation to perfect matching sampling [54], but understanding and controlling the mixing behavior of the Markov chain associated with general Markov bases is much more challenging [25, 61].

5 Conclusion and discussion

In this work, we consider the problem of constructing graphs with prescribed discrete curvatures. We take a first step in addressing this methodological gap in network geometry and discrete curvature by approaching the problem in two ways. First, in Section 3 we develop an evolutionary MCMC-type algorithm to construct graphs whose curvatures approximate a given set. Using this algorithm, we then explore the ensemble of graphs with a fixed discrete curvature sequence. Second, in Section 4 we solve the exact reconstruction problem for Forman–Ricci curvature and degree sequences. We show that there exists a finite set of graph rewiring moves — Markov moves and transpositions — which connect the set of all graphs with a given Forman–Ricci curvature and degree sequence.

Our work opens several new directions for future research. Most results obtained in this article are limited to specific choices of discrete curvatures, but the same questions are equally relevant for all other notions of discrete curvature. More specifically, our algorithmic approach could be extended by optimizing the expensive computations in the mutation/selection phases for specific curvatures. This would make the sampling algorithm practical for a wider range of curvatures and provide a tool to compare among them. Furthermore, a more careful choice of algorithm design might lead to an understanding of and control over the mixing times and convergence rates of the algorithm.

Our results on exact reconstruction for Forman–Ricci target sequences are not as directly extendable to other discrete curvatures. The analysis and solution heavily depends on the simple combinatorial nature of Forman–Ricci curvatures and its connection to joint degrees in a graph. However, an important follow-up question would be to find ways to decrease the prohibitively large size of the Markov bases $\mathcal{M}_{\Delta}$. Adding constraints to the curvature or degree sequences and thus restricting to subspaces of $\ker_{\mathbb{Z}}(A_{\Delta})$ is one possible approach. If the Markov basis size can be controlled, this would enable a practical implementation of the proposed Algorithm 1.
 
Acknowledgements: The authors thank Florentin Münch for a discussion that lead to the proof of Proposition 4.10.

Appendix A Algebraic background for Markov bases

This brief appendix defines the algebraic concepts that appear in the fundamental theorem of Markov bases, Theorem 4.5, and that may not be commonly known. We refer to [17, 25] for further background.

Let $K=\mathbb{C}[x_{1},\dots,x_{n}]$ be the ring of polynomials in variables $x_{1},\dots,x_{n}$ and complex coefficients. For a vector $u=(u_{1},\dots,u_{n})\in\mathbb{N}^{n}$ we write $\mathbf{x}^{u}=x_{1}^{u_{1}}\cdot\dots\cdot x_{n}^{u_{n}}\in K$ for the monomial with exponent vector $u$, and for a vector $\beta\in\mathbb{Z}^{n}$ we write $\beta^{+}$ for the vector with entries $(\beta^{+})_{i}=\beta_{i}$ if $\beta_{i}\geq 0$ and zero otherwise, and $\beta^{-}=\beta^{+}-\beta$.

An ideal in $K$ is a subset $I\subseteq K$ such that $f\in I,g\in K\Rightarrow fg\in I$ and $f,g\in I\Rightarrow f+g\in I$. An ideal $I$ is generated by $f_{1},\dots,f_{k}\subseteq I$ if every $g\in I$ can be written as a polynomial combination $g=\sum_{i=1}^{k}h_{i}f_{i}$ for some $h_{1},\dots,h_{k}\in K$.

The fundamental theorem of Markov bases is thus an algebraic statement which says that Markov bases are precisely the generators of some ideal associated with the kernel $\ker_{\mathbb{Z}}(A)$ of an integer matrix $A\in\mathbb{Z}^{d\times n}$; these are called toric ideals.