2.1 Support size of DERGMs

The number of graphs in the support of a ERGM is $2^{{n \choose 2}}$ . Since a DERGM restricts the support, a natural question that arises is: what is the number of graphs in the support of a DERGM with degeneracy parameter k? Unfortunately, there are no simple formulas to count the number of k-degenerate graphs; nonetheless, we can obtain an asymptotic lower bound as follows.

Theorem 1 (Support size of DERGMs) Let $S_k(n)$ denote the number of simple graphs with n nodes and degeneracy at most k. Then, for a fixed k, there exist positive constants $c_1,c_2 > 0$ and an integer $n_0$ such that for all $n > n_0$ ,

\begin{align*}c_1 \cdot n\;{\log}\;n \leq\;{\log}\;S_k(n) \leq c_2 \cdot n\;{\log}\;n\end{align*}

That is, for a fixed k, and as n goes to infinity, $\log S_k(n) = \Theta\left( n\;{\log}\;n \right).$ On the other hand, for $k=n-1$ , $\log S_{n-1}(n) = \Theta(n^2).$

Theorem 1 is an asymptotic statement that gives an asymptotic upper and lower bound on the support size of DERGMs, when k is a fixed constant. For the finite sample settings, we can consider $k = O(1)$ , that is k is a bounded from above by a constant, whereas n is increasing. (As a practical example, n may be 5,000, but k may be 50 or even 10.) Under such settings, Theorem 1 shows that there are about $O(2^{n\log n})$ graphs in the support of DERGM. On the other hand, the ERGM has $O(2^{n^2})$ graphs. Note that $S_{n-1}(n)$ is the size of the support of the full ERGM. We found two interesting properties: that parameter estimates of a DERGM do not change drastically from that of the corresponding ERGM, see Section 5.2 for a concrete example; and that the graphs eliminated from the support of the ERGM are precisely the ones that cause instability issues, as illustrated in the next result.

Proof of Theorem 1. We derive both upper and lower bounds for the DERGMs support size. A natural lower bound on the number of k-degenerate graphs is the number of well-ordered k-degenerate graphs. A well-ordered k-degenerate graph is a labeled graph with vertex labels $1, \ldots, n$ such that the ordering of the vertices by their labels is a well-ordering of the graph. From Bauer et al. (Reference Bauer, Krug and Wagner2010), the number of well-ordered graphs with degeneracy at most k is given by

\begin{align*} D_k(n) = D_k(n-1) \cdot \sum_{i=0}^{\min\!(n-1,k)} {n-1 \choose i}\end{align*}

By definition, $D_k(n)$ is a lower bound on the $S_k(n)$ . Applying the recursion, for a constant k, we get

\begin{align*}D_k(n) = \left( \sum_{i=0}^{k} {n-1 \choose i} \right) \cdot \left( \sum_{i=0}^{k} {n-2 \choose i} \right) \ldots \cdot \left( \sum_{i=0}^{k} {k \choose i} \right) \cdot \left( \sum_{i=0}^{k-1} {k-1 \choose i} \right) \cdot \left( \sum_{i=0}^{1} {1 \choose i} \right)\end{align*}

which further simplifies as follows:

\begin{align*} D_k (n) &= \prod_{r=k+1}^{n-1} \sum_{i=0}^{k} {r \choose i} \cdot \prod_{r=1}^{k}\sum_{i=0}^{r} {r \choose i}\\ &=\prod_{r=k+1}^{n-1} \sum_{i=0}^{k} {r \choose i} \cdot \prod_{r=1}^{k}2^r\\ &=\prod_{r=k+1}^{n-1} \sum_{i=0}^{k} {r \choose i} \cdot 2^{{k \choose 2}} \end{align*}

Taking logarithms gives

\begin{align*} \;{\log}\;D_k (n) &=\sum_{r=k+1}^{n-1}\;{\log}\;\left(\sum_{i=0}^{k} {r \choose i} \right) + {k \choose 2}\log 2\\ &\geq \sum_{r=k+1}^{n-1}\;{\log}\;{r \choose k} + {k \choose 2}\log 2 \end{align*}

Note that the second term depends only on k, and hence, we can focus on the first term. Let

\begin{align*}T_k(n) \;:\!=\; \sum_{r=k+1}^{n-1}\;{\log}\; {r \choose k} \end{align*}

Using the lower bound ${r \choose k} \geq (r/k)^k$ , we get,

\begin{align*} T_k(n) &\geq k \cdot \sum_{r=k+1}^{n-1}\;{\log}(r/k) \\ &\geq k \cdot \left(\sum_{r=k+1}^{n-1}\;{\log}\;r\right ) - k\;{\log}\;k (n-k-1)\\ &=k \cdot \left(\sum_{r=1}^{n-1}\;{\log}\;r - \sum_{r=1}^k\;{\log}\;r\right ) - k\log k (n-k-1) \\ &= k \cdot \left(\log \!(n-1)! - \;{\log}\;k! \right ) - k\log k (n-k-1) \\ &=\Omega(n\;{\log}\;n) \end{align*}

Thus, the claimed lower bound follows: $\log S_k(n) \geq\;{\log}\;D_k(n) \geq T_k(n) = \Omega(n\;{\log}\;n)$ .

For the upper bound on the support size of k-degenerate graphs, we will use the following strategy. Let $\#G(n,\leq m)$ denote the number of graphs on n nodes with at most m edges, we will show below that

(3) \begin{align}\log \#G(n,\leq m) \leq 2m \cdot {\log}(en)\end{align}

From Proposition 1 below, the maximum number of edges in a k-degenerate graph is $k \cdot$ $n - \left(\substack{(k+1)\\2}\right)$ . Using the fact that

\begin{align*}\mathcal{G}_{n,k} \subset G(n,\leq m)\end{align*}

where $m =k \cdot n - \left(\substack{(k+1)\\2}\right)$ , we have the following upper bound:

\begin{align*}\log S_k(n) &\leq\;{\log}\;\#G\left(n,\leq k \cdot n - {(k+1) \choose 2}\right) \\& \leq 2 \left(k \cdot n - {(k+1) \choose 2}\right)\;{\log}(en) \\& < 2k \cdot n\;{\log}(en) = O(n\;{\log}\;n)\end{align*}

Finally, to see that the upper and lower bounds for the case when $k=n-1$ hold, note that $k=n-1$ is the full ERGM and we have $2^{n \choose 2}$ graphs in the support of an ERGM. Thus $\log S_{n-1}(n) =\;{\log}\;2^{n \choose 2} = \Theta(n^2)$ .

All that remains to be shown is Equation (3). Note that the number of graphs on n nodes with m edges is ${{n \choose 2} \choose m}$ , since there are ${n \choose 2}$ possible locations to choose from and place the m edges. Now the number of graphs with at most m edges is given by

\begin{align*}\#G(n,\leq m) &= \sum_{i=0}^{m} {{n \choose 2} \choose i} \leq \left(\frac{e{n \choose 2}}{m}\right)^m\end{align*}

from the well-known fact $\sum_{i=0}^m {n \choose i} \leq \left(\frac{en}{m}\right)^m$ . Taking logs, we get

\begin{align*}\log \#G(n,\leq m) &\leq\;{\log}\;\left(\frac{e{n \choose 2}}{m}\right)^m \\&\leq\;{\log}\;\left(\frac{en^2}{m}\right)^m \\&\leq 2m\log en\end{align*}

2.2 Stability of sufficient statistics

By restricting the support to include only those graphs with degeneracy at most k, where k is small compared to n, we eliminate “dense” graphs from the model. In turn, this has a stabilizing effect on the sufficient statistics. A formal definition of a stable sufficient statistic in ERGMs is given in Schweinberger (Reference Schweinberger2011) .

Definition 3 (Stable sufficient statistics) Let $S_k(n)$ be the size of support of a DERGM with sufficient statistic t(g). Then, t(g) is said to be stable if for any constant $C > 0$ there exists an integer $n_0$ such that for every $n \geq n_0$

\begin{align*} \underset{g \in \mathcal{G}_{n,k}}{\max} t(g) < C \cdot {\log}\;S_k(n) \end{align*}

or in other words, $\underset{g \in \mathcal{G}_{n,k}}{\max} t(g) \in o(\log S_k(n))$ . On the other hand, t(g) is said to be unstable if for any $C > 0$ , however large,

\begin{align*} \underset{g \in \mathcal{G}_{n,k}}{\max} t(g) \geq C \cdot {\log}\;S_k(n) \end{align*}

A vector of sufficient statistics is stable if all the components of the vector are stable, if any component is unstable, the vector of sufficient statistics is unstable.

Roughly, a sufficient statistic is stable if it can eventually be strictly upper-bounded by the log of the support size of the DERGM. If it cannot be upper-bounded by the log of support size, then it is unstable. For an ERGM, with no support restriction, this definition reduces to strictly upper bounding the sufficient statistic by ${n \choose 2}$ , where n is the number of nodes and it strengthens the definition of stable sufficient statistics in Schweinberger (Reference Schweinberger2011) . The edge-triangle ERGM is not stable due to the instability of the number of triangles, as shown in Schweinberger (Reference Schweinberger2011). However, it turns out that the edge-triangle DERGM is stable.

1. 1. $\underset{g \in \mathcal{G}_{n,k}}{\max} e(g) = k \cdot n - {(k+1) \choose 2}$

2. 2. $\underset{g \in \mathcal{G}_{n,k}}{\max} \triangle(g) = {k \choose 3} + {k \choose 2}(n - k).$

Proof. For this proof, we use the notion of a shell index of a node: define the i-th shell of a graph g to be the difference of the two consecutive cores $H_i(g)\setminus H_{i-1}(g)$ . Note that a node may belong to more than one core, but shell membership is unique. Thus, we say that a vertex v is said to have shell index i if $v\in H_i(g)$ but $v\not\in H_{i+1}(g)$ .

For any given network, the shell sequence $s_1 \leq s_2 \ldots \leq s_n$ is the sorted sequence of shell indices of each node. From Proposition 10 in Karwa et al. (Reference Karwa, Pelsmajer, Petrović, Stasi and Wilburne2017), the maximum number of edges in a graph with a shell sequence $s_1 \leq s_2 \ldots \leq s_n$ is given by:

\begin{align*} {k \choose 2} + \sum_{i=1}^{n-k} s_i\end{align*}

This expression is maximized by graphs in which all the nodes are in the $k{\rm th}$ core, which has a shell sequence $s_1 = k, s_2 = k, \ldots s_n = k$ . Thus, the maximum number of edges in a k-degenerate graph is

\begin{align*} {k \choose 2} + \sum_{i=1}^{n-k} k = \frac{k(k-1)}{2} + k(n-k) = nk - {(k+1) \choose 2}\end{align*}

Similarly, from Proposition 12 in Karwa et al. (Reference Karwa, Pelsmajer, Petrović, Stasi and Wilburne2017), the maximum number of triangles in a graph with shell sequence $s_1 \leq s_2 \ldots \leq s_n$ is given by:

\begin{align*} {k \choose 3} + \sum_{i=1}^{n-k} {s_i \choose 2}\end{align*}

This expression is maximized also when all the nodes are in the $k{\rm th}$ core. Thus, the maximum number of triangles is

\begin{align*} {k \choose 3} + \sum_{i=1}^{n-k} {k \choose 2} = {k \choose 3} + (n-k) {k \choose 2} \\[-38pt] \end{align*}

Proposition 1 shows that the number of triangles in a k-degenerate graph is O(n), whenever $k = O(1) $ . (In fact k can be allowed to grow with n, albeit slowly, see the next theorem) On the other hand, without any restriction on the degeneracy, the number of triangles can be as large as $O(n^3)$ making the ERGMs unstable. The number of triangles in k-degenerate graphs is linear in n, which make them a good candidate to model sparse graphs, which are commonplace in the real world.

In Theorem 2, we use Proposition 1 to show that the edge-triangle DERGM is stable. The way we defined a DERGM assumes that k is fixed; however, note that Theorem 2 shows that k can grow with n, albeit slowly: For instance, if k grows with $\sqrt{\log\!(n)}$ , then the sufficient statistics are still stable.

Proof. We need to show that for all $c>0$ , there exists $n_0$ , there exists $n > n_0$ such that $\max_g (e(g), \triangle(g)) < c \cdot {\log}\;S_k(n))$ where the max is over the support set $g \in \mathcal{G}_{n,k}$ . Fix a g in $\mathcal{G}_{n,k}$ . From Proposition 1, we have,

\begin{align*} (e(g), \triangle(g)) &\leq \left(k \cdot n - {(k+1) \choose 2},{k \choose 3} + {k \choose 2}(n - k)\right) \\ &\leq O(k \cdot n, k^2 \cdot n) \end{align*}

Thus, if $k = o(\sqrt{\log n})$ , we have $(e(g), \triangle(g)) = o(n\log n) = o(\log S_k(n))$ .

2.3 Non-degeneracy of DERGMs

We now show that stability of sufficient statistics implies that a DERGM is non-degenerate. Let us begin by defining degeneracy of a distribution or more precisely the degeneracy of a parameter associated with a distribution. Consider a DERGM defined by the parameter vector $\theta$ and sufficient statistics t(g) and let $M_k(\theta)$ be the set of modes, that is

\begin{align*}M_k(\theta) = \underset{g \in \mathcal{G}_{n,k}}{\arg\max } \frac{e^{\theta^T \cdot t(g)} }{c_k(\theta)}\end{align*}

One also defines a set of $\varepsilon$ -modes for any $0< \varepsilon<1$ :

\begin{align*}M_{\varepsilon,k}(\theta) = \left\{ G \in \mathcal{G}_{n,k} \;:\; e^{\theta^T \cdot t(G)} > (1-\varepsilon) \underset{g \in \mathcal{G}_{n,k}}{\max } e^{\theta^T \cdot t(g)}\right\}\end{align*}

A parameter $\theta$ is said to be asymptotically degenerate if the distribution induced by $\theta$ asymptotically places all of its mass on its modes.

Definition 4 (Asymptotically degenerate parameters, see also Schweinberger Reference Schweinberger2011) A parameter $\theta$ is said to be asymptotically degenerate if

\begin{align*}\lim\limits_{n \rightarrow \infty} \mathbb P_{\theta}(G \in M_k(\theta)) = 1\end{align*}

If, on the other hand, $\lim\limits_{n \rightarrow \infty}\mathbb P_{\theta}(G \in M_k(\theta)) $ is bounded away from 1, the model is asymptotically non-degenerate. We define asymptotic near-degeneracy for DERGMs similarly using $\varepsilon$ -modes.

As Schweinberger (Reference Schweinberger2011) discusses, strict degeneracy in discrete exponential families is not attainable; thus, $\theta$ is said to be near-degenerate if the mass concentrates on $\varepsilon$ -modes. The same reference proves that unstable sufficient statistics lead to near-degenerate distributions. In the following result, we prove that, under a technical condition that the number of graphs in the $\varepsilon$ -modes grows slower than square root of the model support size, stability implies non-(near-)degeneracy in the more general case of DERGMs.

Proof. To show that a DERGM is not near-degenerate, we need to show that $\lim\limits_{n \rightarrow \infty} \mathbb P_{\theta}(G \in M_{\varepsilon,k}(\theta)) < 1$ . That is, we need to show that for every $0<\varepsilon<1$ , however small, $\mathbb P_{\theta}(G \in M_{\varepsilon,k}(\theta)) $ is bounded away from 1 asymptotically.

\begin{align*} \mathbb P_{\theta}(G \in M_{\varepsilon,k}(\theta)) &= \frac{1}{c_k(\theta)} \sum_{g \in M_{k,\varepsilon}(\theta)} \exp\!(\theta^T\cdot t(g)) \\ &= \frac{\sum_{g \in M_{\varepsilon,k}(\theta)}{\exp\!(\theta^T \cdot t(g)) }}{ \sum_{g \in \mathcal{G}_{n,k} }{\exp\!(\theta^T \cdot t(g)) } }\\ &= \frac{\sum_{g \in M_{\varepsilon,k}(\theta)}{\exp\!(\theta^T \cdot t(g)) }}{\sum_{g \in M_{\varepsilon,k}(\theta)} \exp\!(\theta^T \cdot t(g)) + \sum_{g \in \mathcal{G}_{n,k} \setminus M_{\varepsilon,k}(\theta)} \exp\!(\theta^T \cdot t(g))} \\ &=\frac{1}{1 + r_n}\end{align*}

where

\begin{align*} r_n &= \frac{\sum_{g \in \mathcal{G}_{n,k} \setminus M_{\varepsilon,k}(\theta)} e^{\theta^T \cdot t(g)}}{\sum_{g \in M_{\varepsilon,k}(\theta)}{e^{\theta^T \cdot t(g)} }}\end{align*}

Now, showing that $\lim\limits_{n \rightarrow \infty} \mathbb P_{\theta}(G \in M_{\varepsilon,k}(\theta)) < 1$ is equivalent to showing $\lim\limits_{n \rightarrow \infty} r_n > 0$ .

Let $N_m = |M_{\varepsilon,k}(\theta)|$ and let $U_{n,k}(\theta)= \underset{g \in \mathcal{G}_{n,k}}{\max} \theta^T \cdot t(g)$ , and $L_{n,k} = \underset{g \in \mathcal{G}_{n,k}}{\min} \theta^T \cdot t(g)$ . Without loss of generality, we can assume that $L_{n,k}(\theta)$ is 0. This follows from observing that $\mathbb P_{\theta}(G=g)$ is invariant under the translations of $\theta^T \cdot t(g)$ by $-L_{n,k}(\theta)$ . Also, note that for any $g \in M_{\varepsilon,k}(\theta)$ , and any $0 < \varepsilon <~1$ , we have $\theta^T \cdot t(g) \leq U_{n,k}(\theta)$ . Thus, we have

\begin{align*} r_n &= \frac{\sum_{g \in \mathcal{G}_{n,k} \setminus M_k(\theta)} e^{\theta^T \cdot t(g)}}{\sum_{g \in M_{\varepsilon,k}(\theta)}{\exp(\theta^T \cdot t(g)) }} \\ &>\frac{\sum_{g \in \mathcal{G}_{n,k} \setminus M_k(\theta)} e^{\theta^T \cdot t(g)}}{N_m e^{U_{n,k}(\theta)}} \\ & \geq \frac{\sum_{g \in \mathcal{G}_{n,k} \setminus M_k(\theta)}e^{L_{n,k}(\theta)}}{N_m e^{U_{n,k}(\theta)}} = \frac{\sum_{g \in \mathcal{G}_{n,k} \setminus M_k(\theta)}e^{0}}{N_m e^{U_{n,k}(\theta)}} =\frac{S_k(n) - N_m}{N_m e^{U_{n,k}(\theta)}} = \frac{\frac{S_k(n)}{N_m} - 1 }{e^{U_{n,k}(\theta)}} \geq \frac{\frac{S_k(n)}{2N_m} }{e^{U_{n,k}(\theta)}} \\ &\geq \frac{c_0 \sqrt{S_k(n)}}{2e^{U_{n,k}(\theta)}} (\mbox{By assumption, } N_m < c_0 \cdot \sqrt{S_k(n)})\\ &\geq \frac{c_0}{2} \frac{\sqrt{e^{c_1 \cdot n\;{\log}\;n}}}{e^{U_{n,k}(\theta)}} (\mbox{Since }{\log\,}S_k(n) > c_1\cdot n\log n, \mbox{ from Theorem 1}) \end{align*}

The last inequality follows from Theorem 1, which states that there exists a constant $c_1$ , and an $n_0$ such that for all $n > n_0$ , $\log S_k(n) \geq c_1 \cdot n\;{\log}\;n$ . Recall that t(g) being stable means that for all $c> 0$ , there exists an $n_0$ such that for all $n > n_0$ , $\underset{g \in \mathcal{G}_{n,k}}{\max} t(g) < c \cdot {\log}\;S_k(n)$ . Thus, for all $c>0$ ,

\begin{align*}U_{n,k}(\theta) &= \underset{g \in \mathcal{G}_{n,k}}{\max} \theta^T \cdot t(g) \\&< c_{\theta} \cdot c \cdot \;{\log}(S_k(n)) \\&< c_{\theta} \cdot c \cdot c_2 \cdot n\;{\log}\;n\end{align*}

The last inequality again follows from Theorem 1 which states that there exists a constant $c_2$ and an $n_0$ such that for all $n > n_0$ , $\log S_k(n) \leq c_2 \cdot n\;{\log}\;n$ . Here, $c_{\theta}$ is a constant that depends on $\theta$ . Thus, we get, for all $c>0$ , there exists an $n_0$ , $c_1$ and $c_2$ such that for all $n > n_0$ ,

\begin{align*} r_n &> \frac{c_0}{2}\frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{U_{n,k}(\theta)}} \\ &>\frac{c_0}{2} \frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{c \cdot c_{\theta} c_2 \cdot n\;{\log}\;n}}\\ &>\frac{c_0}{2} e^{\left(\frac{c_1}{2} - c \cdot c_2 c_{\theta}\right) \cdot n\;{\log}\;n}\end{align*}

Since this holds for any $c>0$ , let us choose c such that $\frac{c_1}{2} - c\cdot c_2 c_{\theta} = 0$ . Then, $r_n > \frac{c_0}{2}>0$ in the limit, as required.

In order to show an explicit example of a model for which we can find a set of parameter values $\theta$ for which Theorem 3 holds, we spell out the result for the example of the triangle DERGM studied in the previous section. At the same time, we can prove stronger result, relaxing the assumption on the degeneracy k.

1. 1. $k = o(\sqrt{\log n})$ .

For $\theta\in\Theta$ , suppose that:

1. 1. $|\theta|_1 < o(\log n)$ , where $|\theta|_1$ is the $l_1$ norm of $\theta$ ,

2. 2. $\theta\in\Theta$ is such that there exists and constant $c_0$ and an $n_0$ such that for all $n> n_0$ , $|M_{\varepsilon,k}(\theta)| < c_0 \sqrt{S_k(n)}$ , that is the number of graphs in the set of $\varepsilon$ modes does not grow larger than the square root of the total number of graphs in the support of the DERGM.

Then, the edge-triangle DERGM is asymptotically non-(near)-degenerate at $\theta$ .

Assumption 1 of course holds for fixed values of $k;$ thus, it is not restrictive on the DERGM as we defined it, but rather is a relaxation. The last assumption is the same as in the theorem above. Note that the former (concerning the growth of k) is weak, whereas the latter (concerning the number of modes) is strong.

Proof. To prove asymptotic non-(near-)degeneracy, we repeat the same steps as in the theorem above, but consider a finer lower bound on the ratio $r_n$ from the end of the proof:

\begin{align*} r_n > \frac{c_0}{2} \cdot \frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{U_{n,k}(\theta)}} \end{align*}

Now, let us examine $r_n$ for the case of number of edges and triangles. From Proposition 1, there exists a constant $c_2$ and an $n_0$ such that for all $n> n_0$ , the following holds:

\begin{align*} U_{n,k}(\theta) &= \max_g\!(\theta_1, \theta_2)^T \cdot (e(g),\triangle(g)) < |\theta|_1 \cdot \max_g\!(e(g)+\triangle(g)) \\ &< |\theta|_1 \cdot c_2 \cdot k^2 \cdot n \end{align*}

Thus, we have

\begin{align*} r_n &\geq \frac{c_0}{2} \cdot \frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{U_{n,k}(\theta)}} \\ &\geq \frac{c_0}{2} \cdot \frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{|\theta|_1 \cdot c_2\cdot k^2 \cdot n}} \end{align*}

If we allow $|\theta|_1 =o(\;{\log}\;n)$ , and $k = o(\sqrt{\log n})$ , then we have $c_2 |\theta|_1 \cdot k^2 \cdot n = o(n\;{\log}\;n)$ , which means for all $c>0$ , there exists an $n_0$ such that for all $n> n_0$ , $c_2 |\theta|_1 \cdot k^2 \cdot n < c \cdot n\;{\log}\;n$ . Thus, we have

\begin{align*} r_n &\geq \frac{c_0}{2} \cdot \frac{e^{\frac{c_1}{2} \cdot n\;{\log}\;n}}{e^{c \cdot n\;{\log}\;n }} \end{align*}

Choosing $c = \frac{c_1}{2}$ , we get $r_n \geq \frac{c_0}{2}$ , as needed.

Corollary 1 shows that the edge-triangle DERGM is asymptotically non-(near)-degenerate for $k = o(\log n)$ and $|\theta|_1 = o(\log n)$ . This result implies that for large n, the edge-triangle DERGM cannot place all its mass on the set of $\varepsilon$ -modes, and there must be a considerable amount of mass assigned to points outside the set of $\varepsilon$ -modes.

3.1 Sampling graphs from a DERGM with a fixed parameter

In this section, we discuss an MCMC algorithm for sampling graphs from the DERGM for a fixed value of $\theta$ with degeneracy parameter k. The key issue is that to sample from a DERGM using MCMC, we need to ensure that the proposed graphs are in the set $\mathcal G_{n,k}$ , that is they have degeneracy restricted to k. To this end, we consider two different approaches: the first, straightforward approach, is to use the usual tie-no-tie proposal (see, for example, Caimo & Friel Reference Caimo and Friel2011) along with the Metropolis-Hastings step. Such a proposal may generate graphs outside the set $\mathcal G_{n,k}$ , which are naturally rejected by the Metropolis-Hastings algorithm. Thus, whenever the degeneracy of the proposed graph is more than k, the graph is rejected; otherwise, it is accepted with the usual acceptance probability that depends on the change statistics, see Hunter et al. (Reference Hunter, Handcock, Butts, Goodreau and Morris2008a) for more details. Note that the degeneracy of a graph can be computed in O(m) time, where m is the number of edges, using the algorithm of Batagelj & Zaversnik (Reference Batagelj and Zaversnik2003).

While the first method works, it can be wasteful and slow, that is at each step of the Markov chain, we have to compute the degeneracy of the graph and reject it whenever it is larger than k. The second approach is to directly propose graphs from the set $\mathcal G_{n,k}$ . For this, we develop a uniform sampler that proposes graphs uniformly from the set of all k-degenerate graphs. The uniform sampler is presented in Section 6.

Algorithm 1 summarizes the approach 2 where the proposal is the uniform distribution from $\mathcal G_{n,k}$ , denoted by $\mathcal U_{n,k}$ . Let $\pi(g) \propto \exp(\theta_{0}^t t(g))$ . The Metropolis-Hastings acceptance ratio becomes

\begin{align*}\alpha(g_{current}, g_{proposed}) = \min \left(1, \frac{\pi(g_{proposed})}{\pi(g_{current})} \right)\end{align*}

Algorithm 1:

Independent Metropolis algorithm to sample from the model

3.2 Existence of MLE and the approximate MLE

There are two likelihood functions: the true likelihood $l(\theta)$ given by Equation (4) and the estimated likelihood $\hat l(\theta)$ given by Equation (6). Correspondingly, there are two maximizers, the true MLE $\hat \theta$ and the approximate MLE $\tilde \theta$ . We will discuss the existence of the true MLE and the approximate MLE and argue that using a smaller k makes the estimation of the MLE easier.

Using the standard theory of exponential families (Barndorff-Nielsen Reference Barndorff-Nielsen2014), existence of the true MLE $\hat \theta$ depends on the marginal polytope, that is, the convex hull of sufficient statistics of the set $\mathcal{G}_{n,k}$ . The log-likelihood function is concave and a unique maximum exists if and only if the observed sufficient statistic $t(g_{obs})$ lies in the relative interior of the marginal polytope. The marginal polytopes of ERGMs are difficult to obtain in general (see for example Engström & Norén Reference Engström and Norén2011) and known only in few special cases, such as Rinaldo et al. (Reference Rinaldo, Petrović and Fienberg2013), Karwa & SlavkoviĆ (Reference Karwa and Slavković2016). Obtaining the marginal polytopes for the degeneracy-restricted ERGMs appears to be more difficult and is an open problem in general, as it can only be computed for one specific DERGM at a time. We will compute these polytopes numerically for the edge-triangle DERGM in Section 4.

On the other hand, existence of the approximate MLE can be checked numerically. As discussed in Handcock (Reference Handcock2003), the estimated log-likelihood (6) can be written as the log-likelihood of a model from a discrete exponential family with support over $t(G_1), \ldots, t(G_B)$ with observed sufficient statistic $t(g_{obs})$ . Hence, using again the standard theory of exponential families (Barndorff-Nielsen Reference Barndorff-Nielsen2014), one can show that the estimated log-likelihood is concave and Equation (6) has a unique maximum if and only if 0 lies in the interior of the convex hull of $\{\Delta(G_1,g_{obs}), \ldots, \Delta(G_B,g_{obs})\}$ . Thus, assuming that the MLE exists, the existence of the approximate MLE is crucially tied to the sampling algorithm used to approximate the likelihood, which in turn depends on the behavior of the model.

4.1 Insensitivity and lack of degeneracy of DERGMs

We begin by studying the effect of k on the mean value and the natural parameters of DERGMs. The goal is to gain insight into the model degeneracy and excessive sensitivity of DERGMs as a function of k. Roughly, the model is said to suffer from degeneracy issues, if the mean value parameters of the model are pushed to the boundary for different values of the natural parameter. Similarly, the model is said to suffer from excessive sensitivity, small changes in the values of the natural parameters lead to large changes in the mean value parameter, see Schweinberger (Reference Schweinberger2011) for more details.

In the rest of the section, we focus on one-parameter exponential families. We will work with normalized sufficient statistics. Specifically, let $U_k$ denote the maximum of t(g) when $g \in G_{n,k}$ . Let the normalized sufficient statistic be $t_{norm}(g) = t(g)/U_k$ . For the natural parameter $\theta$ , the mean value parameter is given by $ \mu_k(\theta) = \mathbb E_{\mathbb P_{\theta, k}} t_{norm}(g)$ .

We consider two different DERGM models: the two-star DERGM with the number of two-stars as the sufficient statistic, and the triangle DERGM with the number of triangles as the sufficient statistic. Degeneracy corresponds to the situation where if $\theta > 0$ , $\mu_k(\theta) \rightarrow 1$ and $\,\theta <0$ , $\mu_k(\theta) \rightarrow~\!0$ . Sensitivity corresponds to the situation where the derivative of $\mu_k(\theta)$ with respect to $\theta$ is very large in a small neighborhood of $\theta$ .

Recall that our goal is to study the map from $\theta$ to $\mu_k(\theta)$ for varying values of k and gain insights into model behavior. To avoid any issues due to MCMC sampling, we compute this map exactly for a small network, where enumeration is possible. Specifically, we consider networks defined on $n=7$ nodes. When $n=7$ , there are a total of $2^{7 \choose 2}$ possible simple networks. We enumerate all possible networks and compute the number of edges, two-stars, triangles and degeneracy of each network. The total number of networks with different degeneracy values is shown in Table 1.

Table 1.

Number of graphs of degeneracy exactly k for $n=7$ nodes

The plot of mean value vs natural parameter for each DERGM model is generated as follows. We fix a value of k and fix a sufficient statistic. Next, we vary $\theta$ from $-3$ to 3 in steps of $0.01$ . For each value of $\theta$ , we compute the corresponding mean value parameter $\mu_k(\theta)$ using the enumerated networks. We normalize $\mu_k(\theta)$ to make sure it lies between 0 and 1 and plot the normalized $\mu_k(\theta)$ on y-axis and the natural parameter $\theta$ on the x-axis. We repeat this process for different values of k, and obtain a separate plot for each value of k. Similarly, we get different sets of plots for each DERGM. The results are shown in Figures 2 and 3.

Figure 2.

Mean value parameters vs natural parameters for the 2-star DERGM for $n=7$ and $k=2,3,6$ respectively.

Figure 3.

Mean value parameters vs natural parameters for the triangle DERGM for $n=7$ and $k=2,3,6$ respectively.

Let us focus on Figure 3(c). This figure shows the map between $\theta$ and $\mu_k(\theta)$ for the triangle DERGM when $k=6$ and $n=7$ , which is the same as the ERGM (since $k=6$ is the maximum possible, there is no support restriction). The plot shows that the mean value parameter is pushed to its corresponding boundaries for positive and negative values of $\theta$ , that is for $\theta > 0$ , $\mu_k(\theta)$ is close to 1, and for $\theta < 0$ , $\mu_k(\theta)$ is close to 0. Moreover, for $\theta$ close to 0, the mean value parameter is very sensitive to small changes in $\theta$ . This is the classic model degeneracy and excessive sensitivity. On the other hand, if we consider Figure 3(a) and (b), we can see that if we restrict the support to 2-degenerate graphs or 3-degenerate graphs, the mean value map improves. Specifically, for $k=2$ , Figure 3(a) shows that $\mu_k(\theta)$ is not pushed to its boundaries for positive or negative values of $\theta$ , and has a small derivative near $\theta=0$ . This shows that the model does not suffer from degeneracy and excessive sensitivity when k is small. A similar conclusion holds for the 2-star model shown in Figure 2. We also created such plots for $n=50$ , for which we had to resort to MCMC sampling to estimate the mean value parameters, see Figure 4 for the triangle DERGM for $k=3$ and $k=50$ . The results for this setting were the same as described here: For small values of k, the triangle and the two-star DERGM does not suffer from excessive sensitivity and model degeneracy.

Figure 4.

Simulated plot of mean value parameters vs natural parameters, based on MCMC, for the triangle DERGM for $n=50$ and $k=3$ and $k=49$ respectively.

4.2 Existence of approximate MLE, the model polytope, and entropy

Consider first the issue of existence of the approximate MLE. Recall from Section 3.2 that in the MCMC-MLE estimation, the approximate MLE does not exist when the observed sufficient statistics lies outside of the convex hull of the sufficient statistics sampled from $\mathbb P_{\theta_0}$ . In DERGMs, this is more likely to happen when the degeneracy parameter k is large relative to the observed graph degeneracy.

As an example to illustrate this phenomenon, consider fitting the edge-triangle DERGM to Sampson monastery data (Sampson Reference Sampson1968), in particular, the time period T4, available at Batagelj & Mrvar (Reference Batagelj and Mrvar2006) and Hunter et al. (Reference Hunter, Handcock, Butts, Goodreau and Morris2008a). In this dataset, $n=18$ and observed graph degeneracy is $k=3$ . Building on the correspondence between MLE non-existence and the model polytope from Rinaldo et al. (Reference Rinaldo, Fienberg and Zhou2009), we study the location of the observed edge-triangle vector with respect to estimated DERGM model polytopes for varying values of k. Recall that the model polytope for an exponential family model is the convex hull of all observable vectors of sufficient statistics. We estimate the DERGM polytopes as convex hulls of edge-triangle pairs of networks obtained by sampling graphs uniformly from the support $\mathcal{G}_{n,k}$ , using Algorithm 2. Figure 5 shows the estimated model polytopes for different values of k, along with the relative location of the Sampson edge-triangle vector. When $k=3$ , the observed sufficient statistic lies well in the relative interior of the sampled sufficient statistics. On the other hand, when $k=6$ and higher, the observed sufficient statistic lies well outside the convex hull.

Figure 5.

Estimated edge-triangle DERGM model polytopes for increasing k, with x marking the location of the observed sufficient statistic of the Sampson graph. Sample size is 100,000 each; we have verified that the results do not change when sample size is increased to 1,000,000. We do not plot the estimated polytopes for all other values of $k>6$ , but the reader can rest assured that the observed value of the sufficient statistics of the Sampson graph only gets farther removed from the convex hull.

As k increases, the observed edge-triangle count is progressively pushed out of the estimated polytope and becomes probabilistically less likely under the uniform distribution (blue corresponds to lower probability). This is because for larger k, the uniform sampler places more weight on edge-triangle counts of denser graphs, making more sparse edge-triangle counts such as those from Sampson graph probabilistically less likely to appear. Thus, for larger k, the observed edge-triangle count of the Sampson graph lies in the tails of the distribution induced by the uniform sampler. This in turn effects the MCMC-MLE as follows: For larger k, the observed sufficient statistic lies close to the boundary of the true model polytope, or as the figures show, outside the estimated polytope. Unless the MCMC algorithm finds a $\theta_0$ that generates graphs around the observed sufficient statistic, the approximate MLE will not exist. However, this is difficult, since as the observed sufficient statistic approaches the boundary, the number of network configurations corresponding to it becomes smaller. This concept can be formalized by measuring the entropy, discussed next.

Entropy. As explained in Rinaldo et al. (Reference Rinaldo, Fienberg and Zhou2009) (see Section 3.4 therein for details), the shape of the model polytope supports the argument that the full ERGM is ill-behaved. Specifically, they use Shanon’s entropy, which captures the degree to which the model concentrates its mass on network configurations associated with a relatively very small number of network statistics. The rationale is that degenerate models have large areas of low entropy. The correspondence between the model polytope and model degeneracy derived by Rinaldo et al. shows that the extremal rays of the normal fan of the model polytope correspond to directions of the ridges of Shanon’s entropy function where it converges to some fixed value. These extremal rays are outer-normals of the facets (in our case, edges) of the polytope; we see that as k grows, the polytope becomes “flatter” or, equivalently, the directions of the outer-normals of the edges on the lower hull get closer together, making the area of high entropy smaller. Although the exact plots are unavailable for the full ERGM on $n=18$ , we know that for $n=9$ already the rays of normal fan concentrate in a small area of the plane implying that the model has low entropy and is degenerate for a vast majority of parameter values; cf. Rinaldo et al. (Reference Rinaldo, Fienberg and Zhou2009), Figure 4a. As the authors in the said article justify, we use the mean value parameters to illustrate this behavior, where it can be clearly seen.

In contrast, Figure 6 shows that the higher-entropy region is more “spread out” across the parameter region for the DERGMs with smaller values of k. While one cannot, of course, conclude that the model is non-degenerate for all possible parameter values, it is clear that the size of the parameter space that correspond to degenerate regions is certainly less than in the full ERGM. Regarding the caveat that the Figures are also estimated and not exact, we are nevertheless confident in the results, because (1) the algorithm used is a uniform sampler of well-ordered graphs from the model support $\mathcal{G}_{n,k}$ , and (2) the estimated polytope is not far off from the true model polytope: it is missing some extremal graphs that are probabilistically unlikely to be generated by the uniform sampler from the space of graphs $\mathcal{G}_9=\mathcal{G}_{9,8}$ .

Figure 6.

Comparison of degenerate (low entropy) regions in the mean value parameter space for the edge-triangle DERGMs on $n=18$ and increasing k. Sample sizes are 100,000. As k increases, the high-entropy region becomes smaller.

4.3 The likelihood surface changes with k

The shape of the estimated likelihood function changes as we change k. To illustrate this, we use the uniform sampler given in Algorithm 2 to sample graphs uniformly from the support of the full ERGM $\mathcal{G}_n=\mathcal{G}_{n,n-1}$ and $\mathcal{G}_{n,k}$ with $k<n-1$ for various DERGMs, and estimate the likelihood function using the sampled graphs for the Sampson network. Figure 7 shows the contours of the (estimated) likelihood function for various values of $\theta = (\theta_{1},\theta_{2})$ . This figure uncovers an interesting trend: the likelihood surface becomes “flatter” around the maximum value as k grows, making it more difficult to find the maximum itself after a certain number of steps.

Figure 7.

Contour plots of the estimated edge-triangle DERGM likelihood functions for the Sampson network, for various values of $(\theta_1,\theta_2)$ . Here, $n=18$ and $k=3,4,17$ . Note that $k=17$ corresponds to the full ERGM. The estimated likelihood is based on an iid sample of 25,000 graphs in $\mathcal{G}_{n,k}$ .

5.1 Examples where DERGMs fit whereas ERGM fit fails to converge

We first start by showing three examples where the MCMC-MLE procedure fails to converge when fitting an edge-triangle ERGM, whereas it converges when using the edge-triangle DERGM with the degeneracy parameter set to the observed degeneracy. We consider three networks—an undirected version of the Sampson dataset, the Faux Mesa High network and the undirected version of ecoli network, from the ergm package in R. The summary statistics of these networks are given in Table 2. Note that we are not claiming that the edge-triangle DERGM is the best model for these data. Instead, the point is to illustrate that restricting the degeneracy has a direct impact on MCMC-MLE estimation.

Table 2.

Summary of datasets used to fit the edge-triangle DERGMs

The Sampson network has $n=18$ nodes and $m=41$ edges, with an observed degeneracy $k=3$ . The Faux Mesa High network has 205 nodes and 203 edges, and an observed degeneracy of 3. The ecoli network has $n=423$ nodes, $m=519$ edges with a degeneracy $k=3$ . Note that all the networks have a low observed degeneracy. In particular, the ecoli and the faux mesa high networks are very sparse since the degeneracy is very small in comparison to the number of nodes.

While fitting the edge-triangle ERGM to these networks, the MCMC-MLE combined with the step length procedure failed to converge due to model degeneracy; for a detailed study of this model’s degeneracy, see Rinaldo et al. (Reference Rinaldo, Fienberg and Zhou2009). Specifically, the Markov chain started sampling networks whose number of edges and triangles are very far from the observed network, indicating model degeneracy. On the other hand, there were no such issues when fitting the edge-triangle DERGM and the MCMC-MLE combined with the step length procedure converged. The estimated parameter for the edge-triangle DERGMs for these networks is given in Table 3. There are two sources of standard error here, one from the MCMC estimation and another corresponding to the MCMCMLE. The MCMCMLE standard errors are calculated by using an MCMC estimate of the inverse of the estimated fisher information matrix, as described in Hunter & Handcock (Reference Hunter and Handcock2006).

Table 3.

Fitting the edge-triangle DERGM where the edge-triangle ERGM fit fails. The ${}^*$ denotes level of significance, based on the p-values. (The numbers in the parenthesis are the standard errors of the MCMCMLE.)

${}^{*}p<0.001$

5.2 Examples when both ERGM and DERGM fit converges

We now consider cases where the MCMC-MLE procedure is able to fit both an ERGM and a DERGM to the same dataset. In these cases, we show that the parameter estimates obtained from both these models are very close to each other. We fit the edge-triangle DERGMs and ERGM to the florentine dataset. This dataset has $n=16$ vertices and $m=20$ edges, with a degeneracy parameter $k=2$ . We fit DERGMs with increasing values of $k=2,3,\ldots, 15$ . Note that when $k=15$ , the DERGM is equivalent to the edge-triangle ERGM. The parameter estimates are given in Table 4. This table shows that the edge parameter is more or less the same for all the DERGMs and ERGM. The parameter corresponding to the triangles varies, but is within the margin of the standard error.

Table 4.

Fitting DERGM and ERGM to the Florentine data. The ${}^*$ denotes level of significance, based on the p-values. (The numbers in the parenthesis are the standard errors of the MCMCMLE.)

${}^{*}p<0.001$

6.1 A uniform sampler for well-ordered graphs from $\mathcal{G}_{n,k}$

In Bauer et al. (Reference Bauer, Krug and Wagner2010), Algorithm 1, the authors derive a uniform sampler for the set of well-ordered graphs in $\mathcal{G}_{n,k}$ . A well-ordered graph is one in which the node labels are ordered so that no vertex has more than k neighbors with a higher label.

Using this algorithm as a starting point, we make several key changes to ensure that their algorithm is computationally efficient: we convert their algorithm from a recursive one to an iterative one. By doing this, we eliminate many complexity problems inherent in the original algorithm. Specifically, the iterative version eliminates stack overflow issues for large graphs, as well as greatly reduces the execution time of generating a graph.

Let us take a closer look at the following algorithm, based on Bauer et al. (Reference Bauer, Krug and Wagner2010), Algorithm 1, which we improved and updated to a scalable version.

Algorithm 2:

Generate a well-ordered g from $\mathcal{G}_{n,k}$ uniformly.

The algorithm was originally formulated using recursion, which we emulate using two for-loops. The first for-loop populates a list of degrees where each index of the list corresponds to the respective vertex label. The degrees for each vertex are generated using a restricted binomial distribution. Instead of utilizing the cumulative distribution and using binary search to obtain values as suggested by the original paper, we opt to use the probability density function and store the values in a list data structure, reducing the complexity of obtaining the degree values. When the loop reaches the very last vertex, we add that vertex to the working vertex set. For each iteration in the second for-loop, a temporary copy of the current working vertex set is created. We then uniformly generate $d_{i}$ indices to sample without replacement from the vertex set copy and use these samples for the edge set of the current vertex. It is obvious that this sample is uniformly generated, complying with the original algorithm.

For a benchmark, we tested the original recursive version (including generating all possible combinations) and the new iterative version on a machine with the following specifications: Intel Core i7-4790K CPU @ 4.00 GHz, 8 GB DDR3 RAM, Arch Linux x64, with the results shown in Table 5. The results clearly indicate that the scalable version is superior in regard to time complexity. In some applications, it may be desirable to further restrict the sample space of the model by restricting the total number of edges of the graph, or use such a restriction for stratified sampling of $\mathcal{G}_{n,k}$ . To that end, let $\mathcal{G}_{n,m,k}$ be the set of graphs on n nodes and degeneracy k with exactly m edges. Bauer et al. (Reference Bauer, Krug and Wagner2010), Algorithm 2 offer an algorithm for uniform sampling of $\mathcal{G}_{n,m,k}$ ; however, it was not implemented due to the complexity of step 3 that the authors suggest be implemented using Equation (2.7) in Bauer et al. (Reference Bauer, Krug and Wagner2010). Pre-computation of degrees proved nearly impossible in practice for several reasons. The recursive nature of calculating the cardinality for possible graphs of given vertices, edges, and degeneracy yielded very inefficient computations in which the run time of each computation was longer than trying to generate whole graphs by other means. While we were able to alleviate this issue somewhat by utilizing a dynamic programming approach with memoization, even for semi-sparse, average size graphs, numerical overflow occurred, which rendered the speed increase fruitless. Instead, we opt to use Bauer et al. (Reference Bauer, Krug and Wagner2010), Algorithm 3, which is a non-uniform but fast sampler of $\mathcal{G}_{n,m,k}$ . Our implementation of this algorithm, outlined in Algorithm 3, stays true to the pseudo-code given in the original paper, with the only alteration being utilizing the same approach to uniform selection as in our implementation of Algorithm 2.

Table 5.

Run times of the uniform samplers

Algorithm 3:

Generate a well-ordered g from $\mathcal{G}_{n,m,k}$ non-uniformly.

6.2 Stratified sampling of $\mathcal{G}_{n,k}$ to include non-well-ordered graphs if needed

Another issue with Bauer et al. (Reference Bauer, Krug and Wagner2010), Algo.1 is that it generates only so-called “well-ordered” graphs in $\mathcal{G}_{n,k}$ . This misses a part of graphs in the support of our model. To remedy this issue, we classify all missing graphs and produce them via stratified sampling with two strata. Specifically, Algorithm 2 is used to sample from the set of well-ordered graphs in $\mathcal{G}_{n,k}$ , while Algorithm 4, described below, is used to generate non-well-ordered graphs in $\mathcal{G}_{n,k}$ . Let $n_1$ and $n_2$ be the number of well-ordered and non-well-ordered graphs, respectively. The formula for $n_1$ is provided in Bauer et al. (Reference Bauer, Krug and Wagner2010) under the notation $D_n^{(k)}$ , while $n_2$ is studied below. To the best of our knowledge, the literature does not provide a good estimate of the number $n_1$ of well-ordered k-degenerate graphs compared to the total number of k-degenerate graphs. Although we derived a lower bound on the total number of k-degenerate graphs ( $\Omega(n\log n)$ ) in Theorem 1, in this section we study the ratio of $n_1$ and $n_2$ further, which is needed from an algorithmic point of view. It should be noted that, in practice, the uniform sampler from Section 6.1 may only be omitting a tiny fraction of graphs in the support of the DERGM; this situation is described in detail at the end of this Section. Therefore, the reader interested in applications more than in theory behind the algorithms that may not be necessary in practice may skip the remainder of this technical section.

A graph $g\in\mathcal{G}_{n,k}$ is not well-ordered if there exists at least one vertex j with at least $k+1$ neighbors in the set $\{j+1,\dots,n\}$ . Among all such vertices with too many big neighbors, let $k+c$ be the minimum such number of big neighbors, and let i be the index of the smallest vertex that has $k+c$ big neighbors. We construct non-well-ordered graphs and use them to estimate $n_1$ by going through possible cases for the values of c and i. For each case $c=1,\dots,n-k-1$ , some vertex i has $k+c$ neighbors in the set $\{i+1,\dots,n\}$ . For each of the cases, the vertex i can be chosen from the set $\{1,\dots,n-(k+c)\}$ . Note that these $k+c$ neighbors of i can be connected in any arbitrary way, as long as the entire graph is in $\mathcal{G}_{n,k}$ . Thus, we proceed as follows: construct a random graph h on $k+c$ vertices whose labels are in the set $\{i+1,\dots,n\}$ . Then, construct a suspension g over h using vertex i, that is, ensure that i is connected to all $k+c$ vertices of h. Finally, the vertices $\{1,\dots,i\}$ can be connected in any way such that, by minimality of i, the resulting subgraph on $\{1,\dots,i\}$ is well-ordered and, additionally, each vertex in the set $\{1,\dots,i\}$ can have at most k neighbors in the vertex set $\{i+1,\dots,n\}$ . The construction is outlined in Algorithm 4.

Algorithm 4:

Generate a non-well-ordered g from $\mathcal{G}_{n,k}$

There are ${ n-i \choose k+c }$ ways to choose the neighbors of the vertex i on Line 4 and for each choice of neighbors there are $2^{k+c\choose 2}$ graphs $\tilde h$ generated on Line 3. There are $D^{(k)}_i$ well-ordered graphs on Line 7 and $i\sum_{p=1}^k{n-i \choose p}$ graphs on Line 8. Thus, Algorithm 4 constructs the following number of graphs g:

(8) \begin{align} \underbrace{ \sum_{i=1}^{n-(k+1)}\underbrace{{n-i\choose k+1}}_{\mbox{Line 4}} \cdot \underbrace{2^{\left(\substack{k+1\\2}\right)}}_{\mbox{Line 3}} \cdot \underbrace{D_i^{(k)}}_{\mbox{Line 7}} \cdot \underbrace{i\sum_{p=1}^k{n-i \choose p}}_{\mbox{Line 8}} }_{c=1} \nonumber \\ + \underbrace{ \sum_{i=1}^{n-(k+2)}{n-i\choose k+2}\cdot 2^{\left(\substack{k+2\\2}\right)} \cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p} +\cdots }_{c=2} \nonumber \\ \dots + \underbrace{ \sum_{i=1}^{n-(k+n-k-1)}{n-i\choose n-1}\cdot 2^{n-1\choose 2}\cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p}}_{c=n-k-1} \end{align}

(9) \begin{align} = 2^{k+1\choose 2} \cdot\sum_{i=1}^{n-(k+1)}{n-i\choose k+1}\cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p} \nonumber \\ + 2^{k+2\choose 2} \cdot \sum_{i=1}^{n-(k+2)}{n-i\choose k+2}\cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p} +\cdots\nonumber \\ \dots + 2^{n-1\choose 2}\cdot {n-1\choose n-1}\cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p}\end{align}

where each of the $n-k-1$ summands corresponds to one of the cases c.

Note that Equation (9) is an upper bound on $n_2$ , since it counts all graphs g constructed by Algorithm 4. It is also a strict upper bound on the number of graphs g actually returned by the algorithm, since it counts those graphs whose degeneracy happens to be strictly larger than k.

Equation (9) counts all graphs on $k+c$ nodes, $2^{\left(\substack{k+c\\2}\right)}$ , constructed in Step 3. Surely, a better count can be obtained by replacing $2^{\left(\substack{k+c\\2}\right)}$ by

\begin{align*}2^{\left(\substack{k+c\\2}\right)}-\#\{\mbox{well-ordered graphs on $k+c$ vertices of degeneracy}>k\}.\end{align*}

Doing this replacement in the equation is, crucially, still an upper bound on $n_2$ (since the well-ordered graphs of degeneracy larger than k certainly do not contribute to any non-well-ordered graphs of degeneracy at most k). Since

\begin{align*} &\#\{\mbox{well-ordered graphs on $k+c$ nodes of degeneracy}>k\} \\ =&\#\{\mbox{all well-ordered graphs on $k+c$ nodes except those of degeneracy} \leq k\} \\ =& D^{(k+c-1)}_{k+c}-D^{(k)}_{k+c} \end{align*}

the following is a better upper bound on the number of graphs, we wish to keep from Algorithm 4 and thus also an upper bound on $n_2$ :

(10) \begin{align} \sum_{i=1}^{n-(k+1)}{n-i\choose k+1} &\cdot \left( 2^{\left(\substack{k+1\\2}\right)}-\left( D^{(k+1-1)}_{k+1}-D^{(k)}_{k+1}\right)\right) \cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p} + \nonumber \\ \sum_{i=1}^{n-(k+2)}{n-i\choose k+2} &\cdot \left( 2^{\left(\substack{k+2\\2}\right)} -\left( D^{(k+2-1)}_{k+2}-D^{(k)}_{k+2}\right)\right) \cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p} +\cdots \nonumber \\ \dots + {n-1\choose n-1} &\cdot \underbrace{\left( 2^{n-1\choose 2} -\left(D^{(n-1-1)}_{n-1}-D^{(k)}_{n-1}\right)\right) }_{\mbox{Line 3 minus well-ordered of degen}>k}\cdot D_i^{(k)} \cdot i\sum_{p=1}^k{n-i \choose p}\end{align}

Let

\begin{align*}t_{true}=\log n_1/(n_1+n_2)\end{align*}

be the true threshold used to divide the sample in two strata and define

\begin{align*}t_{estimated}=\log n_1/(n_1+(10))\end{align*}

Given that $(10)>n_2$ , $ t_{estimated} < t_{true}\leq 0$ . Therefore, we take the following approach: (1) compute the threshold $t_{estimated}$ for the fixed n and k for which we wish to run the current simulation. (2) If $t_{estimated}$ is close to 0, then that forces $t_{true}$ to be close to 0, which in turn means that there is a very, very small number of non-well-ordered graphs for that choice of n and k and therefore the stratified sampler essentially reduces to sampling well-ordered graphs only.

Of course, if $t_{estimated}$ is not relatively close to 0, then for those values of n and k, while it is possible that $t_{true}$ is close to 0, one should implement both the well-ordered and non-well-ordered algorithm. Falling back on the well-ordered algorithm is equivalent to using an approximate sampler in practice. The users may additionally prefer to replace Algorithm 4 by instead permuting the vertices of the output of Algorithm 2, allowing it to reach the entire sample space $\mathcal{G}_{n,k}$ in another way.