A Python package for statistical validation of weighted bipartite networks using maximum entropy null models.
- Five null models: BiCM, BiWCM, BiECM, BiPECM, BiCReMA
- Binarization utilities: direct thresholding and RCA (Revealed Comparative Advantage)
- Statistical validation with analytical and Monte Carlo p-values
- Multiple testing corrections (Bonferroni, FDR)
- Projection validation for weighted and binary network projections
- Model comparison using information criteria (AIC, AICc, BIC)
- Efficient sampling from null model ensembles
pip install wbnmFor optional dependencies:
pip install "wbnm[viz]" # matplotlib (for plotting)
pip install "wbnm[dev]" # development tools (pytest, black, mypy)
pip install "wbnm[full]" # all dependenciesFor development (editable install from source):
git clone https://github.com/lbuffa/wbnm.git
cd wbnm
pip install -e ".[dev]"import numpy as np
from wbnm import BiCM, BiECM, BiWCM, BiPECM, BiCReMA
from wbnm.validation import fdr
W = np.array([
[5, 0, 3, 1],
[0, 2, 4, 0],
[1, 3, 0, 2],
[0, 1, 2, 0],
[2, 0, 0, 3],
], dtype=float)
model = BiECM(W)
model.solve()
expected = model.expected_matrix()
pvalues = model.pvalues(W)
rejected, adjusted_pvalues = fdr(pvalues, alpha=0.05)
samples = model.sample(n_samples=1000, seed=42)
# GPU is used automatically if available; override with device="cuda" or "cpu"
model_gpu = BiECM(W, device="cuda")All models require every node to have at least one link (non-zero degree/strength). Remove isolated rows/cols before fitting.
Preserves the degree sequence in expectation. Binary model, weights are ignored.
from wbnm import BiCM
from wbnm.utils import rca
import torch
# Direct binarization
B_direct = (W > 0).astype(float)
# RCA binarization (Revealed Comparative Advantage >= 1)
B_rca = rca(torch.as_tensor(W, dtype=torch.float64))
model = BiCM(B_direct)
model.solve()Preserves the strength sequence (sum of weights per node) in expectation.
from wbnm import BiWCM
model = BiWCM(W)
model.solve()Preserves degree sequence and strength sequence in expectation.
from wbnm import BiECM
model = BiECM(W)
model.solve(method="fixed-point", tol=1e-8, max_iter=10000)Preserves the strength sequence and total number of edges in expectation.
from wbnm import BiPECM
model = BiPECM(W)
model.solve()Approximates the BiECM by decoupling topology and weights via two sequential stages:
- Solves a BiCM to reproduce the degree sequence
- Conditioned on those probabilities, solves a BiWCM to reproduce the strength sequence
Same parameter count as BiECM but faster to solve on large networks.
from wbnm import BiCReMA
model = BiCReMA(W)
model.solve(method="fixed-point", beta=0.9)The fixed-point solver supports momentum parameters (beta, beta_schedule, warmup_steps); see the Momentum Parameters table in the API Reference section below.
# Analytical p-values
pvalues = model.pvalues(W)
# Multiple testing correction
from wbnm.validation import fdr, bonferroni
rejected, adjusted = fdr(pvalues, alpha=0.05)from wbnm.validation import pvalues_projection, pvalues_projection_binary
pvals_weighted = pvalues_projection(model, n_samples=1000, rows=True)
pvals_binary = pvalues_projection_binary(model, n_samples=1000, rows=True)
rejected, _ = fdr(pvals_weighted, alpha=0.05)from wbnm.utils import compare_models, print_model_comparison
models = [BiWCM(W).solve(), BiECM(W).solve(), BiPECM(W).solve()]
print_model_comparison(models, criterion="aic")
# criterion: "aic", "aicc", "bic"from wbnm.utils import rca
import torch
# RCA_ij = (w_ij / s_i) / (d_j / W_tot), binary where RCA >= threshold
B = rca(torch.as_tensor(W, dtype=torch.float64), threshold=1.0)from wbnm.utils import biadjacency_to_edgelist, edgelist_to_biadjacency, density
edges = biadjacency_to_edgelist(W, weighted=True)
W_reconstructed = edgelist_to_biadjacency(edges, n_rows=5, n_cols=4)
d = density(W)from wbnm.data import load_example, list_datasets
print(list_datasets())
W, info = load_example("toy")All models share the same interface:
| Method | Description |
|---|---|
solve(method, tol, max_iter, ...) |
Fit the model parameters |
sample(n_samples, seed) |
Generate random samples from the null model |
expected_matrix() |
Expected weight or connection-probability matrix |
log_likelihood() |
Log-likelihood of the observed data under the model |
pvalues(observed) |
Analytical p-values for each entry |
aic(), aicc(), bic() |
Information criteria |
All models expose these properties (computed from the input matrix):
| Property | Description |
|---|---|
n_rows, n_cols |
Matrix dimensions |
degree_sequence_rows/cols |
Number of links per node |
strength_sequence_rows/cols |
Sum of weights per node |
is_solved |
Whether the model has been fitted |
All models support three solvers via the method parameter of solve():
method= |
Description |
|---|---|
"fixed-point" |
Fixed-point iteration (default, stable) |
"quasi-newton" |
Diagonal Newton step, faster near the solution |
"coordinate" |
Coordinate-wise bisection (~1e-18 precision) |
Common parameters: tol (default 1e-8), max_iter (default 10000), verbose, alpha (step size).
The chunk parameter processes the matrix in row-blocks to reduce peak memory on GPU:
model.solve(method="fixed-point", chunk=256)The following parameters tune the fixed-point update rule. They are active only for
method="fixed-point" on BiECM, BiPECM, and BiCReMA. Passing non-default values on
BiCM or BiWCM, or with any other method, raises a UserWarning.
| Parameter | Default | Supported by |
|---|---|---|
beta |
0.9 |
BiECM, BiPECM, BiCReMA |
beta_schedule |
None ("warmup" or "linear") |
BiECM, BiPECM, BiCReMA |
warmup_steps |
500 |
BiECM, BiPECM, BiCReMA |
use_momentum |
True |
BiPECM only |
Models follow the device of the input:
| Input | Device used |
|---|---|
| numpy array / list | CPU |
| torch tensor | same device as tensor |
explicit device="cuda" |
CUDA (any input type) |
pytest tests/ -v
pytest tests/ --cov=wbnm --cov-report=html- Python >= 3.9
- PyTorch >= 1.10.0
- NumPy >= 1.21.0
- Saracco, F., et al. (2015). Randomizing bipartite networks: the case of the World Trade Web. Scientific Reports, 5, 10595.
- Bruno, M., et al. (2023). Inferring comparative advantage via entropy maximization. Journal of Physics: Complexity, 4, 045011.
- Squartini, T., & Garlaschelli, D. (2011). Analytical maximum-likelihood method to detect patterns in real networks. New Journal of Physics, 13, 083001.
- Vallarano, N., et al. (2021). Fast and scalable likelihood maximization for Exponential Random Graph Models with local constraints. Scientific Reports, 11, 15227.
MIT License, see LICENSE for details.
- Lorenzo Buffa - lorenzo.buffa@cref.it
- Daniele Cirulli - daniele.cirulli@cref.it
Centro Ricerche Enrico Fermi (CREF), Rome, Italy