ergmx

Exponential-family random graph models (ERGMs) in Python, with a Rust core.

ergmx fits, simulates, summarizes and checks ERGMs with a high-level API in the spirit of R's ergm and statnet: R-style formulas, the same term names and statistics, and summary() and gof() that read like R's.

74 terms and 9 operators for directed, undirected and bipartite networks; curved ERGMs; sample space constraints; missing ties; multilevel networks (as MPNet); samples of networks (as ergm.multi); temporal ERGMs and dynamic simulation (as tergm); MPLE, contrastive divergence and Monte Carlo MLE; MCMC diagnostics, log-likelihoods, model comparison and goodness of fit; tie probabilities, marginal effects and tables of results, all validated against R. See what's missing.

import ergmx
from ergmx import datasets

g = datasets.load("faux.mesa.high")   # an igraph.Graph; ergm's networks are bundled
fit = ergmx.ergm(
    g, "edges + nodefactor('Sex') + nodematch('Grade') + nodematch('Race') + gwesp(0.5, fixed=TRUE)",
    seed=1,
)
fit.summary()

Monte Carlo Maximum Likelihood Results:

                   Estimate  Std. Error  MCMC %  z value  Pr(>|z|)
edges               -6.1846      0.1734       0  -35.669    <1e-04 ***
nodefactor.Sex.M    -0.1256      0.0747       0   -1.681   0.09274 .
nodematch.Grade      1.9710      0.1758       0   11.210    <1e-04 ***
nodematch.Race       0.2657      0.1188       0    2.235   0.02539 *
gwesp.fixed.0.5      1.2166      0.0853       0   14.271    <1e-04 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Log-likelihood: -867.1531 (MC SE 0.194)   AIC: 1744.3062   BIC: 1784.0461
Converged after 6 iterations (4 chains, 1024 samples).

R's ergm gives -6.1884, -0.1293, 1.9761, 0.2699, 1.2178 and a log-likelihood of -867.56 on the same model, in 16.8 s; ergmx took 2.5 s. (A high-precision estimate of the log-likelihood is -867.09: R's is off by 0.47, see validation.)

Then check the fit, as with R's gof():

result = fit.gof()   # degree, edgewise shared partners, geodesic distances, model statistics
print(result["degree"])
result.plot()

Check the MCMC and compare models, as with R's mcmc.diagnostics() and anova():

fit.mcmc_diagnostics()        # effective sizes, R-hat across chains, Geweke; .plot() for traces
simpler = ergmx.ergm(g, "edges + nodefactor('Sex') + nodematch('Grade') + nodematch('Race')")
ergmx.compare(simpler, fit)   # log-likelihoods, AIC, BIC, likelihood-ratio test

Features

• Networks: igraph.Graph or networkx.Graph/DiGraph, directed or undirected. Vertex attributes are available to the terms; edges with na=True mark missing dyads.

• Formulas in R syntax ("edges + gwesp(0.5, fixed=TRUE)", even with net ~ in front), or terms combined with +: edges() + gwesp(0.5, fixed=True). Nothing is evaluated: arguments must be literals.

• Terms, with ergm's definitions and names:

  	Undirected	Directed
  Dyadic	edges, edgecov, sociality	edges, edgecov, mutual, asymmetric, sender, receiver
  Degree	kstar(k), degree(d), isolates, concurrent, twopath, gwdegree	istar(k), ostar(k), idegree(d), odegree(d), isolates, twopath, gwidegree, gwodegree
  Triads and cycles	triangle, cycle(k), gwesp, gwdsp, gwnsp, esp(d), dsp(d), nsp(d)	triangle, ttriple, ctriple, transitive, cycle(k), and the shared partner terms with any type (OTP, ITP, RTP, OSP, ISP)
  Attributes	nodematch (diff=TRUE too), nodemix, nodefactor, nodecov, absdiff, absdiffcat	the same, plus nodeifactor, nodeofactor, nodeicov, nodeocov
  Operators	offset(term) (with -inf to forbid ties), F(~terms, ~filter)	the same

  Bipartite networks have b1star, b1degree, gwb1degree, b1concurrent, b1factor, b1cov, b1nodematch, b1dsp, gwb1dsp and their b2 twins. The geometrically weighted terms take a fixed decay (fixed=TRUE) or, as in ergm by default, estimate it: curved ERGMs. edgecov('name') reads an n x n matrix from a graph attribute.

• Constraints, as in ergm: bd (bounded degrees, for fixed-choice designs), blocks (fix the dyads of some mixing types), degrees, odegrees and idegrees, with degree-preserving MCMC moves.

• Missing ties: likelihood inference conditional on the observed dyads (Handcock and Gile 2010), assuming they are missing at random, for the estimates, standard errors, log-likelihood and goodness of fit.

• Bipartite networks (bipartite=True): only the ties between modes are modeled, with their own terms and ergm's goodness of fit statistics.

• Curved ERGMs: the decays of gwesp, gwdegree and the other geometrically weighted terms estimated in the MPLE, contrastive divergence, the Monte Carlo MLE and the log-likelihood, with an overflow statistic where ergm's cutoff stops with an error.

• Multilevel networks as one network with a level attribute, as MPNet models them (Wang et al. 2013): ergm's S() for terms within a level or between two (S(~edges + gwesp(0.5, fixed=TRUE), ~level == 'A')), MPNet's cross-level configurations (star2ax, txax, atxax, l3axb, c4axb...), and nodemix, F() and blocks to model kinds of ties.

• Samples of networks, as R's ergm.multi: ergmx.Networks(g1, g2, ...) models many networks (classrooms, households) jointly, and N(~terms, lm=~log(n) + weekday) lets the coefficients depend on network-level attributes, with R's syntax and names.

• Temporal ERGMs, as R's tergm: ergmx.tergm([wave1, wave2, wave3], "Form(~edges + mutual) + Persist(~edges)") fits the conditional MLE of a series of networks, with Form(), Persist(), Diss(), Cross() and Change(); fit.simulate(time_slices=20) and ergmx.simulate_dynamic() run the process forward, with the ties that form and dissolve and their durations.

• Estimation:

  • dyad-independent models: the exact MLE (logistic regression), with log-likelihood, AIC and BIC;
  • other models: Monte Carlo MLE starting from the MPLE (or from the contrastive divergence estimate, init="CD"), with Hummel et al. (2012) step lengths, the log-normal approximation, an adaptive MCMC interval that targets an effective sample size, and standard errors that include the MCMC error;
  • estimate="MPLE" or estimate="CD" for those estimates only;
  • degenerate models stop with a DegeneracyError that says why: a density guard (as in ergm) and a check that the estimate is still moving, with the simulated and observed statistics side by side.

• Log-likelihood of dyad-dependent models, with its Monte Carlo standard error, by path sampling from the dyad-independent submodel as ergm does, but integrated with the Euler-Maclaurin corrected trapezoidal rule: its error falls as 1/bridges^4 instead of 1/bridges^2 for ergm's midpoint rule. ergmx.compare(fit1, fit2, ...) tabulates log-likelihoods, AIC, BIC and likelihood-ratio tests of nested models.

• MCMC diagnostics: fit.mcmc_diagnostics() reports mean deviations from the observed statistics, naive and time-series standard errors, effective sizes, split R-hat across the parallel chains and Geweke z-scores, with trace and density plots.

• MCMC in Rust: tie/no-tie (TNT) proposals mixed with triadic proposals, which close or open triangles, for models with triangle or shared partner terms, directed or not (like ergm's MH_SPDyad default). Chains run in parallel threads.

• Goodness of fit: fit.gof() or ergmx.gof(network, formula, coef) compares the degree (in- and out-degree if directed), edgewise shared partner and geodesic distance distributions, and the model statistics, with simulated networks: tables with Monte Carlo p-values like R's, and plot().

• ergmx.simulate(network, formula, coef, nsim) returns graphs of the same kind as the input, or their statistics; fit.simulate() uses the estimates.

• Interpreting and reporting: fit.predict() (tie probabilities, as ergm's predict()), fit.marginal_effects() (as ergMargins), fit.odds_ratios(), fit.confint(), and ergmx.table(fit1, fit2), a table of models identical to texreg's screenreg(), also as LaTeX, HTML and Markdown; to_frame() for pandas.

• ergmx.summary_stats(network, formula): R's summary(net ~ formula).

• ergmx.datasets: the networks of R's ergm documentation (flomarriage, flobusiness, samplk1-3, faux.mesa.high, faux.dixon.high, faux.magnolia.high), ergm.multi's 318 household networks Goeyvaerts, and multinets' multilevel linked_sim.

Documentation

A user guide, a term reference, the API reference and a guide for R users, built with Sphinx in docs/; every example runs when it is built:

uv sync --group docs
uv run sphinx-build -W --keep-going -d docs/_build/doctrees docs docs/_build/html

.github/workflows/docs.yml publishes it to GitHub Pages.

Validation against R

scripts/r_reference.R fits the models below with ergm 4.12 (ergm.multi 0.3.0 and tergm 4.2.2 for samples and series of networks) and stores the results in tests/data/r_reference.json; the test suite compares.

The networks are ergm's flomarriage, samplk3, faux.mesa.high and faux.dixon.high (248 students, directed friendship nominations), the second and third also with missing dyads, and multinets' linked_sim. R's own standard errors vary by about 15% between seeds on the small networks, so the SE comparison is only as tight as R allows. Full details: docs/validation.md.

Performance

benchmarks/benchmark.R and benchmarks/benchmark.py: median of 3 seeds on an Apple M4 Pro, both with their defaults, which include the log-likelihood. R's ergm runs on one thread; the single-threaded ergmx run is limited to one thread too.

ergmx's log-likelihood samples about 17 times more than ergm's, which is what makes it accurate. Without the log-likelihood on either side (eval_loglik=False and eval.loglik = FALSE), a single thread was 1.9 to 3.6 times faster than R on these models.

The proposal matters as much as the language: with plain TNT proposals, R takes 128 s on faux.magnolia.high.

Architecture

python/ergmx/         Python API
  formula.py            R-style formula parsing (Python's ast, no eval)
  terms.py              term definitions: names, parameters, directedness
  _estimation.py        MPLE, contrastive divergence, Monte Carlo MLE, degeneracy checks
  _loglik.py            log-likelihood by path sampling
  _diagnostics.py       MCMC diagnostics
  _gof.py               goodness of fit
  _compare.py           model comparison
  constraints.py        sample space constraints
  datasets.py           ergm's and multinets' networks, bundled in data/
  _fit.py               ErgmFit and its summary table
  _simulate.py          ergm(), simulate(), summary_stats()
src/                  Rust core (PyO3), exposed as ergmx._core.Model
  network.rs            sorted neighbour lists + edge list for O(1) random ties
  terms.rs              the Term trait and the change statistics
  sampler.rs            Metropolis-Hastings: TNT, triadic and degree-preserving moves
  space.rs              the sample space: free dyads, degree bounds
  lib.rs                bindings; chains run in parallel with rayon

Adding a term means implementing its change statistic in src/terms.rs and describing it in python/ergmx/terms.py.

Not yet

• bd() bounds by attribute, and more of ergm's constraints and terms.
• Valued networks (ergm.count) and egocentric data (ergm.ego).
• tergm's EGMME estimator and duration terms; ergm.multi's gofN(), and N()'s subset, offset and label arguments.

Installation

ergmx needs Python 3.11 or newer. With uv:

uv add "ergmx[igraph,plot]"        # in a uv project; or: uv pip install / pip install

Releases come as wheels with the Rust core already compiled, for Linux, macOS and Windows, so no Rust compiler is needed. To install it from GitHub, which builds it from source:

uv add "ergmx[igraph,plot] @ git+https://github.com/neylsoncrepalde/ergmx"

Building from source needs a C linker (Xcode Command Line Tools, build-essential or the Visual Studio C++ Build Tools) and Rust 1.88 or newer. If Rust isn't installed, maturin downloads a private copy (about 500 MB, cached). The installation guide has the details.

Development

With uv and rustup:

git clone https://github.com/neylsoncrepalde/ergmx.git && cd ergmx
rustup toolchain install    # the Rust pinned in rust-toolchain.toml
uv sync                     # Python 3.14, locked dependencies, ergmx built in release mode
uv run pytest               # Python tests
cargo clippy --release --all-targets -- -D warnings && cargo test --release

uv run rebuilds the Rust core after changes to it. uv build makes a wheel and a source distribution in dist/; .github/workflows/release.yml builds the wheels of every platform and publishes them to PyPI when a GitHub release is published.

To regenerate the R reference results or rerun the benchmark (needs R with ergm, igraph and jsonlite):

Rscript scripts/r_reference.R && Rscript scripts/r_gof_reference.R
Rscript benchmarks/benchmark.R && uv run python benchmarks/benchmark.py

License

GPL-3.0, like R's ergm.