WBNM - Weighted Bipartite Null Models

A Python package for statistical validation of weighted bipartite networks using maximum entropy null models.

Features

• 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

Installation

git clone https://github.com/lbuffa/wbnm.git
cd wbnm
pip install -e .

pip install -e ".[dev]"   # development tools
pip install -e ".[viz]"   # matplotlib
pip install -e ".[full]"  # all dependencies

Quick Start

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")

Null Models

All models require every node to have at least one link (non-zero degree/strength). Remove isolated rows/cols before fitting.

BiCM - Bipartite Configuration Model

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()

BiWCM - Bipartite Weighted Configuration Model

Preserves the strength sequence (sum of weights per node) in expectation.

from wbnm import BiWCM

model = BiWCM(W)
model.solve()

BiECM - Bipartite Enhanced Configuration Model

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)

BiPECM - Bipartite Partial Enhanced Configuration Model

Preserves the strength sequence and total number of edges in expectation.

from wbnm import BiPECM

model = BiPECM(W)
model.solve()

BiCReMA - Bipartite Conditional Reconstruction Method A

Approximates the BiECM by decoupling topology and weights via two sequential stages:

1. Solves a BiCM to reproduce the degree sequence
2. 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.

Validation

# Analytical p-values
pvalues = model.pvalues(W)

# Multiple testing correction
from wbnm.validation import fdr, bonferroni
rejected, adjusted = fdr(pvalues, alpha=0.05)

Projection Validation

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)

Model Comparison

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"

Utilities

Binarization

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)

Network Conversions

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)

Example Datasets

from wbnm.data import load_example, list_datasets

print(list_datasets())
W, info = load_example("toy")

API Reference

Methods

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

Properties

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

Solver Methods

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)

Momentum Parameters (fixed-point only)

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

Device Handling

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)

Testing

pytest tests/ -v
pytest tests/ --cov=wbnm --cov-report=html

Requirements

• Python >= 3.9
• PyTorch >= 1.10.0
• NumPy >= 1.21.0

References

• 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.

License

MIT License, see LICENSE for details.

Authors

• Lorenzo Buffa - lorenzo.buffa@cref.it
• Daniele Cirulli - daniele.cirulli@cref.it

Centro Ricerche Enrico Fermi (CREF), Rome, Italy