Simulates data from a Gaussian Graphical Model (GGM) with a small-world
network structure. The generative process proceeds in three stages. First,
a ring lattice is constructed with neighbors + 1 nearest-neighbor
connections per node, where neighbors is derived from nodes
and density. Second, the lattice is randomly pruned to the target
density, introducing degree heterogeneity from the outset. Third,
edges are rewired with probability rewire, where rewired edges are
placed preferentially on node pairs with higher combined degree (degree-
weighted rewiring). This approach produces more realistic degree
distributions than standard Watts-Strogatz rewiring while preserving
the local clustering structure of the lattice. Edge weights are assigned
by structural priority: edges retained from the pruned lattice receive
larger partial correlation weights based on their original neighbor
distance, grounded in the empirical observation that shorter-distance
connections tend to carry stronger weights in psychometric networks.
The resulting network is used to generate multivariate normal (or skewed)
data. Parameters do not have default values (except
negative_proportion, snr, target_condition,
max_correlation, and max_iterations) and must each be set.
See Details and Examples to get started.
simulate_smallworld(
nodes,
density,
rewire,
snr = 1,
negative_proportion,
sample_size,
skew = 0,
skew_range = NULL,
target_condition = 30,
max_correlation = 0.8,
max_iterations = 100
)A named list with four elements:
Numeric matrix of dimension sample_size x nodes
containing the simulated observations drawn from the population GGM.
Rows are cases; columns are variables named V01, V02,
etc. Values are continuous (or skewed continuous when skew != 0).
To produce ordinal data, pass the columns through
categorize.
A list of input, derived, and estimated parameters:
nodes --- Number of nodes (as supplied)
density --- Target edge density (as supplied)
neighbors --- Number of nearest neighbors k derived
from nodes and density; the initial ring lattice uses
neighbors + 1 connections before pruning
rewire --- Edge rewiring probability (as supplied)
negative_proportion --- Proportion of negative edges,
either as supplied or as sampled from the empirical distribution
sample_size --- Number of simulated observations
skew --- Named numeric vector of per-variable skew values
actually applied (after rounding and possible resampling)
weibull --- Named numeric vector of length 2 giving the
Weibull shape and scale parameters of the absolute edge
weight distribution actually used. If ridge conditioning was applied,
these are re-estimated from the conditioned network via MLE.
omega --- Smallworldness omega statistic of the generated
network (Telesford et al., 2011). Values near zero indicate small-world
structure; negative values indicate lattice-like structure; positive
values indicate random-like structure (see smallworldness)
Q --- Newman-Girvan modularity of the population network
(Omega) computed via igraph::modularity on absolute edge
weights, using the community partition returned by
igraph::cluster_leiden with objective_function =
"modularity". Because cluster_leiden is a heuristic algorithm,
Q reflects a high-quality but not necessarily optimal partition.
Values near zero indicate absence of community structure, consistent
with the network theory of psychopathology.
Population-level network parameters:
R --- Population correlation matrix derived from the GGM
via pcor2cor
Omega --- Population partial correlation matrix (the GGM
edge weight matrix), with zeros for absent edges
Iteration and conditioning diagnostics:
iterations --- Number of sampling attempts needed to find
a valid network
rejections --- Character vector recording the rejection
reason for each failed attempt. Common reasons include disconnected
graph structure, condition number exceeding target_condition,
and maximum correlation exceeding max_correlation
lambda --- Ridge regularization parameter \(\lambda\)
added to the diagonal of the precision matrix to ensure positive
definiteness; NA if no conditioning was required
condition --- Condition number of the final population
correlation matrix R
Numeric (length = 1).
Number of nodes in the network.
Minimum of three nodes. The total number of nodes should be between 8
and 54 to remain within the range of the empirical networks used to fit
the Weibull parameter model; values outside this range are accepted but
will trigger extrapolation and may produce a warning from
weibull_parameters.
Numeric (length = 1).
Target edge density of the network after pruning the initial ring lattice.
Controls the number of nearest neighbors k each node is connected
to via \(k = \mathrm{round}((\text{nodes} \times
\frac{\text{nodes}-1}{2} \times \text{density}) / \text{nodes})\), subject
to a minimum of 2 and a maximum of
\(\lfloor (\text{nodes}-1)/2 \rfloor\). A ring lattice with
neighbors + 1 connections per node is first generated and then
pruned to this density, introducing degree heterogeneity before rewiring.
Must be between 0 and 1. A minimum density sufficient to maintain a
connected graph is enforced; values below this threshold will produce an
informative error.
Numeric (length = 1). Probability of rewiring each edge. Unlike the standard Watts-Strogatz model, rewired edges are placed using degree-weighted random selection: node pairs with higher combined degree have a greater probability of receiving a rewired edge (see Details). Values near 0 preserve the pruned lattice structure; values near 1 produce approximately random networks. The small-world regime typically occurs at intermediate values (roughly 0.01 to 0.30). Must be between 0 and 1.
Numeric (length = 1).
Signal-to-noise ratio of the absolute partial correlation weights,
defined as \(\bar{|w|} / \mathrm{SD}(|w|)\). Values less than 1
produce wider, more heterogeneous weight distributions; values greater
than 1 produce narrower, more homogeneous distributions. Empirically
observed SNR values ranged from 0.648 to 1.712; values outside this
range are accepted but will trigger a warning. Defaults to 1.
Numeric (length = 1). Proportion of edges that are negative (inhibitory). Must be between 0 and 0.50. The upper bound of 0.50 is a mathematical constraint of the sign-flipping procedure used to assign negative edges (see Details). If not provided, a value is sampled from a truncated normal distribution reflecting the empirical distribution of true negative partial correlations across 194 psychometric networks: mean = 0.34, SD = 0.086, bounded to \([0.083, 0.50]\).
Numeric (length = 1).
Number of observations to generate from the population multivariate
normal distribution. Also influences the predicted Weibull scale
parameter via weibull_parameters: larger samples are
associated with smaller, more precisely estimated edge weights.
Numeric (length = 1 or nodes).
Skew applied to each variable after generation from the multivariate
normal. Can be a single value (applied to all variables) or one value
per variable. Values are rounded to the nearest increment of 0.05 in
the range \([-2, 2]\). Defaults to 0 (no skew).
Numeric (length = 2).
If provided, a skew value is drawn uniformly from this range for each
variable, overriding skew. Both values must be between -2 and 2.
Numeric (length = 1).
Target condition number (using kappa with exact =
TRUE) applied when ridge regularization is needed to recover a positive
definite precision matrix. The smallest ridge penalty \(\lambda\) that
brings the condition number to this target is found via root-finding
(uniroot), subject to a maximum shrinkage of approximately
23% following Peeters et al. (2020). Lower values produce
better-conditioned matrices. Defaults to 30. Values up to
100 are accepted but not recommended.
Numeric (length = 1).
Maximum allowed absolute pairwise correlation in the population
correlation matrix R. Any draw where
max(abs(R[lower.tri(R)])) > max_correlation is rejected and a
new attempt is made. Must be between 0 and 1. Defaults to 0.80.
Numeric (length = 1).
Maximum number of attempts to find (1) a connected network structure,
(2) a network satisfying the small-world screening criterion, and (3) a
valid set of edge weights. The error message reports a frequency table
of rejection reasons to assist with diagnosing convergence failures.
Defaults to 100.
Alexander P. Christensen <alexpaulchristensen@gmail.com>
Lattice generation and pruning
The generative process begins by constructing a ring lattice with
neighbors + 1 nearest-neighbor connections per node via
sample_smallworld with p = 0. The
+1 overshoot ensures the lattice always has more edges than the
target density, guaranteeing that pruning can proceed. The lattice is
then randomly pruned to the target density by removing edges
uniformly at random, subject to the constraint that the pruned graph
remains connected. Because removal is random, different nodes lose
different numbers of edges, producing a heterogeneous degree distribution
before any rewiring occurs. This heterogeneity provides a non-uniform
prior for the subsequent degree-weighted rewiring step.
Degree-weighted rewiring
Each edge in the pruned lattice is independently selected for rewiring
with probability rewire. For each selected edge \((i, j)\),
node \(i\) is kept fixed and the \(j\) endpoint is redirected to a
new target node. Valid targets are restricted to node pairs involving
node \(i\) that (1) are not currently connected, and (2) have never
been occupied during the current rewiring pass (i.e., were absent in
the original pruned lattice). Among valid targets, the new endpoint is
selected with probability proportional to \(\sqrt{d_i + d_k}\), where
\(d_i\) and \(d_k\) are the current degrees of the two nodes in
the candidate pair. The square-root transformation moderates the
rich-get-richer tendency of linear preferential attachment, producing
degree heterogeneity consistent with the empirical range of psychometric
networks without generating extreme hubs. Node degrees are updated
incrementally after each rewire so that subsequent rewiring steps
reflect the current graph state.
Edge weight assignment
Edge weights are assigned by ranking edges according to their distance
in the pruned lattice before rewiring. Edges retained from the pruned
lattice receive their original neighbor distance (1 = nearest neighbor,
2 = second nearest, etc.); rewired edges are assigned a distance of
neighbors + 1, placing them at the bottom of the priority order.
Absolute edge weights are drawn from a Weibull distribution whose
parameters are predicted from the network size, sample size, and signal-
to-noise ratio using a Seemingly Unrelated Regression (SUR) model fitted
to 194 empirical psychometric networks (Huth et al., 2025). The sorted
(descending) Weibull weights are then mapped onto the distance ranking
so that the largest weights go to the shortest-distance edges. This
assignment is grounded in the empirical observation that shorter-distance
local connections tend to carry larger partial correlation weights in
psychometric networks.
Negative edge assignment
Negative edges are introduced by flipping the signs of a subset of
nodes — multiplying all edges incident to those nodes by -1. The number
of nodes to flip is derived from negative_proportion via the
inversion formula \((1 - \sqrt{1 - 2p}) / 2\), where \(p\) is the
target proportion of negative edges. This formula assumes that an edge
is negative if and only if exactly one of its endpoints is flipped,
giving an expected negative proportion of \(2 \cdot (k/n) \cdot
(1 - k/n)\) for \(k\) flipped nodes out of \(n\) total. The formula
requires \(p \leq 0.50\), which is why negative_proportion is
bounded at 0.50.
Seminal introduction to the Watts-Strogatz small-world model
Watts, D. J., & Strogatz, S. H. (1998).
Collective dynamics of 'small-world' networks.
Nature, 393(6684), 440--442.
Logic for weight assignments
Muldoon, S. F., Bridgeford, E. W., & Bassett, D. S. (2016).
Small-world propensity and weighted brain networks.
Scientific Reports, 6(1), 22057.
Analytical approximation of random-graph average path length
Newman, M. E. J., Strogatz, S. H., & Watts, D. J. (2001).
Random graphs with arbitrary degree distributions and their applications.
Physical Review E, 64(2), 026118.
Empirical network data used to fit the Weibull SUR model
Huth, K. B. S., Haslbeck, J. M. B., Keetelaar, S., Van Holst, R. J., &
Marsman, M. (2025). Statistical evidence in psychological networks.
Nature Human Behaviour.
Maximum ridge shrinkage bound
Peeters, C. F., van de Wiel, M. A., & van Wieringen, W. N. (2020).
The spectral condition number plot for regularization parameter evaluation.
Computational Statistics, 35(2), 629--646.
# Basic small-world network (moderate density, moderate rewiring)
result <- simulate_smallworld(
nodes = 20,
density = 0.30,
rewire = 0.20,
sample_size = 500
)
# Lattice-like structure (low rewiring preserves local connectivity)
result <- simulate_smallworld(
nodes = 20,
density = 0.30,
rewire = 0.01,
sample_size = 500
)
# Random-like structure (high rewiring destroys lattice regularity)
result <- simulate_smallworld(
nodes = 20,
density = 0.30,
rewire = 0.80,
sample_size = 500
)
# Fix the proportion of negative edges
result <- simulate_smallworld(
nodes = 20,
density = 0.30,
rewire = 0.20,
sample_size = 500,
negative_proportion = 0.20
)
Run the code above in your browser using DataLab