LVBlockState#

class graph_tool.inference.LVBlockState(s, g=None, self_loops=True, sigma=1, **kwargs)[source]#

Bases: DynamicsBlockStateBase

Base state for network reconstruction based on dynamical or behavioral data, using the stochastic block model as a prior. This class is not meant to be instantiated directly, only indirectly via one of its subclasses. Nevertheless, its paramteres are inherited, and are documented as follows.

Parameters:
sndarray of shape (N,M) or list of VertexPropertyMap or VertexPropertyMap

Time series or independent samples used for reconstruction.

If the type is ndarray, it should correspond to a (N,M) data matrix with M samples for all N nodes.

If the parameter g is provided, this can be optionally a list of of VertexPropertyMap objects, where each entry in this list must be a VertexPropertyMap with type vector<int> or vector<double>, depending on wether the model has discrete or continuous state values. If a single property map is given, then a single time series is assumed.

If the parameter t below is given, each property map value for a given node should contain only the states for the same points in time given by that parameter.

gGraph (optional, default: None)

Initial graph state. If not provided, an empty graph will be assumed.

tlist of VertexPropertyMap or VertexPropertyMap (optional, default: [])

If nonempty, this allows for a compressed representation of the time-series parameter s, corresponding only to points in time where the state of each node changes. Each entry in this list must be a VertexPropertyMap with type vector<int> containing the points in time where the state of each node change. The corresponding state of the nodes at these times are given by parameter s. If a single property map is given, then a single time series is assumed.

xEdgePropertyMap or float (optional, default: 1.)

Initial value of the edge weights. If a float is given, the edge weights will be assumed to be the same for all edges.

x_rangepair of float``s (optional, default: ``(-np.inf, np.inf))

Determines the allowed range of edge weights.

thetaEdgePropertyMap or float (optional, default: 0.)

Initial value of the node parameters. If a float is given, the node weights will be assumed to be the same for all edges.

theta_rangepair of float``s (optional, default: ``(-np.inf, np.inf))

Determines the allowed range of the node parameters.

disable_xdistboolean (optional, default: False)

If True the, quantization of the edge weights will be turned off, and \(L_1\) regularization will be used instead.

disable_tdistboolean (optional, default: False)

If True the, quantization of the node parameters will be turned off, and \(L_1\) regularization will be used instead.

nestedboolean (optional, default: True)

If True, a NestedBlockState will be used, otherwise BlockState.

state_argsdict (optional, default: {})

Arguments to be passed to NestedBlockState or BlockState.

bstateNestedBlockState or BlockState (optional, default: None)

If passed, this will be used to initialize the block state directly.

self_loopsbool (optional, default: False)

If True, it is assumed that the inferred graph can contain self-loops.

state_argsint (optional, default: 1 << 16)

Maximum latent edge multiplicity allowed.

nmax_extendint (optional, default: 8)

Maximum number of range expansions during bisection search.

entropy_args: ``dict`` (optional, default: ``{}``)

Override default arguments for entropy() method and releated operations.

Methods

add_edge(u, v, x[, dm])

Add edge \((u, v)\) with multiplicity dm and weight x.

bisect_t(v[, maxiter, tol, entropy_args, ...])

Perform a bisection search to find the best value of node v.

bisect_x(u, v[, maxiter, tol, entropy_args, ...])

Perform a bisection search to find the best weight value for edge \((u, v)\).

collect_marginal([g])

Collect marginal inferred network during MCMC runs.

collect_marginal_multigraph([g])

Collect marginal latent multigraph during MCMC runs.

copy(**kwargs)

Return a copy of the state.

edge_MI(u, v)

Return the mutual information between nodes \(u\) and \(v\), according to their time-series.

edge_TE(u, v)

Return the transfer entropy between nodes \(u\) and \(v\), according to their time-series.

edge_cov(u, v[, toffset, pearson])

Return the covariance (or Pearson correlation if pearson == True) between nodes \(u\) and \(v\), according to their time-series.

edge_mcmc_sweep([beta, niter, k, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample latent edges.

edge_multiflip_mcmc_sweep([beta, niter, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection merge-split MCMC to sample discrete edge weight categories.

entropy([latent_edges, density, aE, sbm, ...])

Return the description length, i.e. the negative joint log-likelihood.

get_block_state()

Return the underlying block state, which can be either BlockState or NestedBlockState.

get_candidate_edges([k, r, max_rk, epsilon, ...])

Return the \(\lfloor\kappa N\rceil\) best edge candidates according to a stochastic second neighbor search.

get_dyn_state([s])

Return an LVState instance corresponding to the inferred model, optionally with initial state given by s.

get_edge_prob(u, v, x[, entropy_args, epsilon])

Return conditional posterior log-probability of edge \((u,v)\).

get_edges_prob(elist[, entropy_args, epsilon])

Return conditional posterior log-probability of an edge list, with shape \((E,2)\).

get_entropy_args()

Return the current default values for the parameters of the function entropy(), together with other operations that depend on them.

get_graph()

Return the current inferred graph.

get_params(params)

Gets the model parameters via the dictionary params.

get_theta()

Return latent node values.

get_tvals()

Return latent node categories.

get_x()

Return latent edge weights.

get_xvals()

Return latent edge weight categories.

mcmc_sweep([beta, niter, edge, edge_swap, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample latent edges and network partitions.

quantize_x(delta)

Quantize weight values according to multiples of \(\Delta\).

remove_edge(u, v[, dm])

Remove edge \((u, v)\) with multiplicity dm.

reset_entropy_args()

Reset the current default values for the parameters of the function entropy(), together with other operations that depend on them.

sample_delta(minval, maxval[, beta, ...])

Sample from the conditional posterior of the \(\Delta\) parameter associated with the edge and node categories.

sample_t(v[, beta, maxiter, tol, ...])

Sample a value for node v according to the conditional posterior.

sample_tl1(minval, maxval[, beta, maxiter, ...])

Sample from the conditional posterior of the \(\lambda\) parameter associated with the node categories.

sample_val_lprob(x, xc[, beta])

Compute probability of sampling value x from bisection history xc and inverse temperature beta.

sample_x(u, v[, beta, maxiter, tol, ...])

Sample a weight value for edge \((u, v)\) according to the conditional posterior.

sample_xl1(minval, maxval[, beta, maxiter, ...])

Sample from the conditional posterior of the \(\lambda\) parameter associated with the edge categories.

sbm_mcmc_sweep([multiflip])

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample node partitions.

set_params(params)

Sets the model parameters via the dictionary params.

set_state(g, w)

Set all edge multiplicities via EdgePropertyMap w.

swap_mcmc_sweep([beta, niter, k, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to swap edge endpoints.

theta_mcmc_sweep([beta, niter, pold, pnew, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample node parameters.

theta_multiflip_mcmc_sweep([beta, pmerge, ...])

Perform sweeps of a Metropolis-Hastings acceptance-rejection merge-split MCMC to sample discrete node value categories.

tvals_sweep([beta, niter, maxiter, ...])

Perform sweeps of a greedy update on the node category values, based on bisection search.

update_edge(u, v, nx)

update edge \((u, v)\) with weight nx.

update_entropy_args(**kwargs)

Update the default values for the parameters of the function entropy() from the keyword arguments, in a stateful way, together with other operations that depend on them.

update_node(v, nt)

update node \((u, v)\) with value nt.

virtual_add_edge(u, v, x[, dm, entropy_args])

Return the difference in description length if edge \((u, v)\) would be added with multiplicity dm and weight x.

virtual_remove_edge(u, v[, dm, entropy_args])

Return the difference in description length if edge \((u, v)\) with multiplicity dm would be removed.

virtual_update_edge(u, v, nx[, entropy_args])

Return the difference in description length if edge \((u, v)\) would take a new weight nx.

virtual_update_node(v, nt[, entropy_args])

Return the difference in description length if node v would take a new value nt.

xvals_sweep([beta, niter, maxiter, tol, ...])

Perform sweeps of a greedy update on the edge weight category values, based on bisection search.

add_edge(u, v, x, dm=1)#

Add edge \((u, v)\) with multiplicity dm and weight x.

bisect_t(v, maxiter=0, tol=1e-07, entropy_args={}, reversible=False, fb=False)#

Perform a bisection search to find the best value of node v.

Parameters:
uint or Vertex

Source

vint or Vertex

Target

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

reversibleboolean (optional, default: False)

Perform search in a manner that is usable for a reversible Makov chain.

fbboolean (optional, default: False)

Perform a Fibonacci (a.k.a. golden ratio) search, instead of a random bisection search.

bisect_x(u, v, maxiter=0, tol=1e-07, entropy_args={}, reversible=False, fb=False)#

Perform a bisection search to find the best weight value for edge \((u, v)\).

Parameters:
uint or Vertex

Source

vint or Vertex

Target

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

reversibleboolean (optional, default: False)

Perform search in a manner that is usable for a reversible Makov chain.

fbboolean (optional, default: False)

Perform a Fibonacci (a.k.a. golden ratio) search, instead of a random bisection search.

collect_marginal(g=None)#

Collect marginal inferred network during MCMC runs.

Parameters:
gGraph (optional, default: None)

Previous marginal graph.

Returns:
gGraph

New marginal graph, with internal edge EdgePropertyMap "eprob", containing the marginal probabilities for each edge.

Notes

The posterior marginal probability of an edge \((i,j)\) is defined as

\[\pi_{ij} = \sum_{\boldsymbol A}A_{ij}P(\boldsymbol A|\boldsymbol D)\]

where \(P(\boldsymbol A|\boldsymbol D)\) is the posterior probability given the data.

collect_marginal_multigraph(g=None)#

Collect marginal latent multigraph during MCMC runs.

Parameters:
gGraph (optional, default: None)

Previous marginal multigraph.

Returns:
gGraph

New marginal graph, with internal edge EdgePropertyMap "w" and "wcount", containing the edge multiplicities and their respective counts.

Notes

The mean posterior marginal multiplicity distribution of a multi-edge \((i,j)\) is defined as

\[\pi_{ij}(w) = \sum_{\boldsymbol G}\delta_{w,G_{ij}}P(\boldsymbol G|\boldsymbol D)\]

where \(P(\boldsymbol G|\boldsymbol D)\) is the posterior probability of a multigraph \(\boldsymbol G\) given the data.

copy(**kwargs)#

Return a copy of the state.

edge_MI(u, v)#

Return the mutual information between nodes \(u\) and \(v\), according to their time-series.

edge_TE(u, v)#

Return the transfer entropy between nodes \(u\) and \(v\), according to their time-series.

edge_cov(u, v, toffset=True, pearson=False)#

Return the covariance (or Pearson correlation if pearson == True) between nodes \(u\) and \(v\), according to their time-series.

edge_mcmc_sweep(beta=inf, niter=1, k=1, elist_args={}, keep_elist=False, pold=1, pnew=1, pxu=0.1, pm=1, premove=1, maxiter=0, tol=1e-07, binary=True, deterministic=False, sequential=True, parallel=True, verbose=False, entropy_args={}, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample latent edges.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

kint (optional, default: 1)

\(\kappa\) parameter to be passed to get_candidate_edges().

elist_argsdict (optional, default: {})

Paramters to pass to call get_candidate_edges().

keep_elistboolean (optional, default: False)

If True, the candidate edge list from last call will be re-used (if it exists).

poldfloat (optional, default: 1)

Relative probability of proposing a new edge weight from existing categories.

pnewfloat (optional, default: 1)

Relative probability of proposing a new edge weight from a new categories.

pxufloat (optional, default: .1)

Probability of choosing from an existing category uniformly at random (instead of doing a bisection search).

pmfloat (optional, default: 1)

Relative probability of doing edge multiplicity updates.

premovefloat (optional, default: 1)

Relative probability of removing edges.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

binaryboolean (optional, default: True)

If True, the latent graph will be assumed to be a simple graph, otherwise a multigraph.

deterministicboolean (optional, default: False)

If True, the the order of edge updates will be determinisitc, otherwise uniformly at random.

sequentialboolean (optional, default: True)

If True, a sweep will visit every edge candidate once, otherwise individiual updates will be chosen at random.

parallelboolean (optional, default: True)

If True, the updates are performed in parallel, using locks on edges candidate incident on the same node.

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

edge_multiflip_mcmc_sweep(beta=inf, niter=1, pmerge=1, psplit=1, pmergesplit=1, pmovelabel=1, gibbs_sweeps=1, c=0.1, maxiter=0, tol=1e-07, accept_stats=None, verbose=False, entropy_args={}, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection merge-split MCMC to sample discrete edge weight categories.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

pmergefloat (optional, default: 1)

Relative probability of merging two discrete categories.

psplitfloat (optional, default: 1)

Relative probability of splitting two discrete categories.

pmergesplitfloat (optional, default: 1)

Relative probability of simultaneoulsly merging and splitting two discrete categories.

pmovelabelfloat (optional, default: 1)

Relative probability of moving the value of a discrete category.

gibbs_sweepsint (optional, default: 1)

Number of Gibbs sweeps performed to achieve a split proposal.

cdouble (optional, default: .1)

Probability of choosing a category uniformly at random to perform a merge, otherwise an adjacent one is chosen.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

accept_statsdict (optional, default: None)

If provided, the dictionary will be updated with acceptance statistics.

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

entropy(latent_edges=True, density=False, aE=1, sbm=True, xdist=True, tdist=True, xl1=1, tl1=1, alpha=1, delta=1e-08, normal=False, mu=0, sigma=1, **kwargs)#

Return the description length, i.e. the negative joint log-likelihood.

Warning

The default arguments of this function are overriden by those obtained from get_entropy_args(). To update the defaults in a stateful way, update_entropy_args() should be called.

Parameters:
latent_edgesboolean (optional, default: True)

If True, the adjacency term of the description length will be included.

densityboolean (optional, default: False)

If True, a geometric prior for the total number of edges will be included.

aEdouble (optional, default: 1)

If density=True, this will correspond to the expected number of edges according to the geometric prior.

sbmboolean (optional, default: True)

If True, SBM description length will be included.

xdistboolean (optional, default: True)

If True, the quantized edge weight distribution description length will be included.

tdistboolean (optional, default: True)

If True, the quantized node parameter distribution description length will be included.

xl1float (optional, default: 1)

Specifies the \(\lambda\) parameter for \(L_1\) regularization for the edge weights if xdist == False, or the Laplace hyperprior for the discrete categories if xdist == True.

tl1float (optional, default: 1)

Specifies the \(\lambda\) parameter for \(L_1\) regularization for the node paraemters if tdist == False, or the Laplace hyperprior for the discrete categories if tdist == True.

deltafloat (optional, default: 1e-8)

Specifies the precision parameter for the discrete categories.

normalboolean (optional, default: False)

If True, a normal distribution will be used for the weight priors.

mudouble (optional, default: 0)

If normal == True, this will be the mean of the normal distribution.

sigmadouble (optional, default: 1)

If normal == True, this will be the standard deviation of the normal distribution.

Notes

The “entropy” of the state is the negative log-likelihood of the generative model for the data \(\boldsymbol S\), that includes the inferred weighted adjacency matrix \(\boldsymbol{X}\), the node parameters \(\boldsymbol{\theta}\), and the SBM node partition \(\boldsymbol{b},\) given by

\[\begin{split}\begin{aligned} \Sigma(\boldsymbol{S},\boldsymbol{X},\boldsymbol{\theta}|\lambda_x,\lambda_{\theta},\Delta) = &- \ln P(\boldsymbol{S}|\boldsymbol{X},\boldsymbol{\theta})\\ &- \ln P(\boldsymbol{X}|\boldsymbol{A},\lambda_x, \Delta)\\ &- \ln P(\boldsymbol{A},\boldsymbol{b})\\ &- \ln P(\boldsymbol{\theta}, \lambda_{\theta}, \Delta). \end{aligned}\end{split}\]

The term \(P(\boldsymbol{S}|\boldsymbol{X},\boldsymbol{\theta})\) is given by the particular generative model being used and \(P(\boldsymbol{A},\boldsymbol{b})\) by the SBM. The weight ditribution is given by the quantized model

\[P(\boldsymbol X|\boldsymbol A,\lambda_x,\Delta) = \frac{\prod_{k}m_{k}!\times \mathrm{e}^{-\lambda_x \sum_k |z_k|}(\mathrm{e}^{\lambda\Delta} - 1)^{K}} {E!{E-1 \choose K-1}2^{K}\max(E,1)}\]

where \(\boldsymbol z\) are the \(K\) discrete weight categories, and analogously

\[P(\boldsymbol\theta|\lambda_{\theta},\Delta) =\frac{\prod_{k}n_{k}!\times \mathrm{e}^{-\lambda \sum_k |u_k|} \sinh(\lambda_{\theta}\Delta)^{K_{\theta}-\mathbb{1}_{0\in\boldsymbol u}} (1-\mathrm{e}^{-\lambda_{\theta}\Delta})^{\mathbb{1}_{0\in\boldsymbol u}}} {N!{N-1 \choose K_{\theta}-1}N},\]

is the node parameter quantized distribution. For more details see [peixoto-network-2024].

References

[peixoto-network-2024]

Tiago P. Peixoto, “Network reconstruction via the minimum description length principle”, arXiv: 2405.01015

[peixoto-scalable-2024]

Tiago P. Peixoto, “Scalable network reconstruction in subquadratic time”, arXiv: 2401.01404

get_block_state()#

Return the underlying block state, which can be either BlockState or NestedBlockState.

get_candidate_edges(k=1, r=1, max_rk='k', epsilon=0.01, c_stop=False, max_iter=0, knn=False, gradient=None, h=1e-06, f_max_iter=10, tol=1e-06, allow_edges=False, include_edges=True, use_hint=True, nrandom=0, keep_all=False, exact=False, return_graph=False, keep_iter=False, entropy_args={}, verbose=False)#

Return the \(\lfloor\kappa N\rceil\) best edge candidates according to a stochastic second neighbor search.

Parameters:
kfloat (optional, default: 1)

\(\kappa\) parameter.

rfloat (optional, default: 1)

Fraction of second neighbors to consider during the search.

max_rkfloat (optional, default: "k")

Maximum number of second-neighbors to be considered per iteration. A string value "k" means that this will match the number of first neighbors.

epsilonfloat (optional, default: .01)

Convergence criterion.

c_stopboolean (optional, default: False)

If True, the clustering coefficient will be used for the convergence criterion.

max_iterint (optional, default: 0)

Maximum number of iterations allowed (0 means unlimited).

knnboolean (optional, default: False)

If True, the KNN graph will be returned.

gradientboolean (optional, default: None)

Whether to use the gradient to rank edges. If None, it defaults to True is the number of edge categories is empty.

hfloat (optional, default: 1e-8)

Step length used to compute the gradient with central finite difference.

allow_edgesboolean (optional, default: False)

Permit currently present edges to be included in the search.

use_hintboolean (optional, default: True)

Use current edges as a hint during the search.

nrandomint (optional, default: 0)

Add this many random entries to the list.

keep_allboolean (optional, default: False)

Keep all entries seen during the search, not only the best.

exactboolean (optional, default: False)

If True an exact quadratic algorithm will be used.

return_graphboolean (optional, default: False)

If True the result will be returned as graph and a property map.

keep_iterboolean (optional, default: False)

If True the result contain the iteration at which an entry has been found.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
elist:class:~numpy.ndarray of shape (E, 2)

Best entries.

a:class:~numpy.ndarray

Edge scores.

get_dyn_state(s=None)[source]#

Return an LVState instance corresponding to the inferred model, optionally with initial state given by s.

get_edge_prob(u, v, x, entropy_args={}, epsilon=1e-08)#

Return conditional posterior log-probability of edge \((u,v)\).

get_edges_prob(elist, entropy_args={}, epsilon=1e-08)#

Return conditional posterior log-probability of an edge list, with shape \((E,2)\).

get_entropy_args()#

Return the current default values for the parameters of the function entropy(), together with other operations that depend on them.

get_graph()#

Return the current inferred graph.

get_params(params)#

Gets the model parameters via the dictionary params.

get_theta()#

Return latent node values.

get_tvals()#

Return latent node categories.

get_x()#

Return latent edge weights.

get_xvals()#

Return latent edge weight categories.

mcmc_sweep(beta=inf, niter=1, edge=True, edge_swap=True, edge_multiflip=True, theta=True, theta_multiflip=True, sbm=True, xvals=True, tvals=True, k=1, keep_elist=False, verbose=False, elist_args={}, edge_mcmc_args={}, edge_swap_mcmc_args={}, edge_multiflip_mcmc_args={}, xvals_mcmc_args={}, theta_mcmc_args={}, theta_multiflip_mcmc_args={}, tvals_mcmc_args={}, sbm_mcmc_args={}, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample latent edges and network partitions.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

kint (optional, default: 1)

\(\kappa\) parameter to be passed to get_candidate_edges().

elist_argsdict (optiona, default: {})

Paramters to pass to call get_candidate_edges().

keep_elistboolean (optional, default: False)

If True, the candidate edge list from last call will be re-used (if it exists).

edgeboolean (optiona, default: True)

Whether to call edge_mcmc_sweep().

edge_mcmc_argsdict (optiona, default: {})

Paramters to pass to call edge_mcmc_sweep().

edge_swapboolean (optiona, default: True)

Whether to call swap_mcmc_sweep().

edge_mcmc_argsdict (optiona, default: {})

Paramters to pass to call swap_mcmc_sweep().

edge_multiflipboolean (optiona, default: True)

Whether to call edge_multiflip_mcmc_sweep().

edge_multiflip_mcmc_argsdict (optiona, default: {})

Paramters to pass to call edge_multiflip_mcmc_sweep().

thetaboolean (optiona, default: True)

Whether to call theta_mcmc_sweep().

theta_mcmc_argsdict (optiona, default: {})

Paramters to pass to call theta_mcmc_sweep().

theta_multiflipboolean (optiona, default: True)

Whether to call theta_multiflip_mcmc_sweep().

theta_multiflip_mcmc_argsdict (optiona, default: {})

Paramters to pass to call theta_multiflip_mcmc_sweep().

sbmboolean (optiona, default: True)

Whether to call sbm_mcmc_sweep().

sbm_mcmc_argsdict (optiona, default: {})

Paramters to pass to call sbm_mcmc_sweep().

xvalsboolean (optiona, default: True)

Whether to call xvals_sweep().

xvals_mcmc_argsdict (optiona, default: {})

Paramters to pass to call xvals_sweep().

tvalsboolean (optiona, default: True)

Whether to call tvals_sweep().

tvals_mcmc_argsdict (optiona, default: {})

Paramters to pass to call tvals_sweep().

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

**kwargsdict (optional, default: {})

Remaining keyword parameters will be passed to all individual MCMC functions.

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

quantize_x(delta)#

Quantize weight values according to multiples of \(\Delta\).

remove_edge(u, v, dm=1)#

Remove edge \((u, v)\) with multiplicity dm.

reset_entropy_args()#

Reset the current default values for the parameters of the function entropy(), together with other operations that depend on them.

sample_delta(minval, maxval, beta=inf, maxiter=0, tol=1e-07, entropy_args={})#

Sample from the conditional posterior of the \(\Delta\) parameter associated with the edge and node categories.

Parameters:
minvalfloat

Minimum value to consider.

minvalfloat

Maximum value to consider.

betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of iterations.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
deltafloat

Sampled value of \(\Delta\).

sample_t(v, beta=inf, maxiter=0, tol=1e-07, entropy_args={}, fb=False)#

Sample a value for node v according to the conditional posterior.

Parameters:
uint or Vertex

Source

vint or Vertex

Target

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

reversibleboolean (optional, default: False)

Perform search in a manner that is usable for a reversible Makov chain.

fbboolean (optional, default: False)

Perform a Fibonacci (a.k.a. golden ratio) search, instead of a random bisection search.

sample_tl1(minval, maxval, beta=inf, maxiter=0, tol=1e-07, entropy_args={})#

Sample from the conditional posterior of the \(\lambda\) parameter associated with the node categories.

Parameters:
minvalfloat

Minimum value to consider.

minvalfloat

Maximum value to consider.

betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of iterations.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
tl1float

Sampled value of \(\lambda\).

sample_val_lprob(x, xc, beta=inf)#

Compute probability of sampling value x from bisection history xc and inverse temperature beta.

sample_x(u, v, beta=inf, maxiter=0, tol=1e-07, entropy_args={}, fb=False)#

Sample a weight value for edge \((u, v)\) according to the conditional posterior.

Parameters:
uint or Vertex

Source

vint or Vertex

Target

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

reversibleboolean (optional, default: False)

Perform search in a manner that is usable for a reversible Makov chain.

fbboolean (optional, default: False)

Perform a Fibonacci (a.k.a. golden ratio) search, instead of a random bisection search.

sample_xl1(minval, maxval, beta=inf, maxiter=0, tol=1e-07, entropy_args={})#

Sample from the conditional posterior of the \(\lambda\) parameter associated with the edge categories.

Parameters:
minvalfloat

Minimum value to consider.

minvalfloat

Maximum value to consider.

betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of iterations.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
xl1float

Sampled value of \(\lambda\).

sbm_mcmc_sweep(multiflip=True, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample node partitions. The remaining keyword parameters will be passed to mcmc_sweep() or multiflip_mcmc_sweep(), if multiflip=True.

set_params(params)#

Sets the model parameters via the dictionary params.

set_state(g, w)#

Set all edge multiplicities via EdgePropertyMap w.

swap_mcmc_sweep(beta=inf, niter=1, k=1, elist_args={}, keep_elist=False, pmove=1, ptmove=1, pswap=1, deterministic=False, sequential=True, parallel=True, verbose=False, entropy_args={}, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to swap edge endpoints.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

kint (optional, default: 1)

\(\kappa\) parameter to be passed to get_candidate_edges().

elist_argsdict (optional, default: {})

Paramters to pass to call get_candidate_edges().

keep_elistboolean (optional, default: False)

If True, the candidate edge list from last call will be re-used (if it exists).

pmovefloat (optional, default: 1)

Relative probability of swaping the weights between two randomly chosen edges.

ptmovefloat (optional, default: 1)

Relative probability of moving a single edge endpoint of an edge with a candidate edge.

pswapfloat (optional, default: 1)

Relative probability of swapping the endpoints of two randomly selected edges.

deterministicboolean (optional, default: False)

If True, the the order of edge updates will be determinisitc, otherwise uniformly at random.

sequentialboolean (optional, default: True)

If True, a sweep will visit every edge candidate once, otherwise individiual updates will be chosen at random.

parallelboolean (optional, default: True)

If True, the updates are performed in parallel, using locks on edges candidate incident on the same node.

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

theta_mcmc_sweep(beta=inf, niter=1, pold=1, pnew=1, maxiter=0, tol=1e-07, deterministic=False, sequential=True, parallel=True, verbose=False, entropy_args={}, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection sampling MCMC to sample node parameters.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

poldfloat (optional, default: 1)

Relative probability of proposing a new node value from existing categories.

pnewfloat (optional, default: 1)

Relative probability of proposing a new node value from a new categories.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

deterministicboolean (optional, default: False)

If True, the the order of node updates will be determinisitc, otherwise uniformly at random.

sequentialboolean (optional, default: True)

If True, a sweep will visit every node once, otherwise individiual updates will be chosen at random.

parallelboolean (optional, default: True)

If True, the updates are performed in parallel.

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

theta_multiflip_mcmc_sweep(beta=inf, pmerge=1, psplit=1, pmergesplit=1, pmovelabel=0, gibbs_sweeps=1, c=0.1, niter=1, maxiter=0, tol=1e-07, entropy_args={}, accept_stats=None, verbose=False, **kwargs)#

Perform sweeps of a Metropolis-Hastings acceptance-rejection merge-split MCMC to sample discrete node value categories.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 1)

Number of sweeps.

pmergefloat (optional, default: 1)

Relative probability of merging two discrete categories.

psplitfloat (optional, default: 1)

Relative probability of splitting two discrete categories.

pmergesplitfloat (optional, default: 1)

Relative probability of simultaneoulsly merging and splitting two discrete categories.

pmovelabelfloat (optional, default: 1)

Relative probability of moving the value of a discrete category.

gibbs_sweepsint (optional, default: 1)

Number of Gibbs sweeps performed to achieve a split proposal.

cdouble (optional, default: .1)

Probability of choosing a category uniformly at random to perform a merge, otherwise an adjacent one is chosen.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

accept_statsdict (optional, default: None)

If provided, the dictionary will be updated with acceptance statistics.

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

tvals_sweep(beta=inf, niter=100, maxiter=0, min_size=1, tol=1e-07, entropy_args={})#

Perform sweeps of a greedy update on the node category values, based on bisection search.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 100)

Number of categories to update.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

min_sizeint (optional, default: 1)

Minimum size of node categories that will be updated.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.

update_edge(u, v, nx)#

update edge \((u, v)\) with weight nx.

update_entropy_args(**kwargs)#

Update the default values for the parameters of the function entropy() from the keyword arguments, in a stateful way, together with other operations that depend on them.

Values updated in this manner are preserved by the copying or pickling of the state.

update_node(v, nt)#

update node \((u, v)\) with value nt.

virtual_add_edge(u, v, x, dm=1, entropy_args={})#

Return the difference in description length if edge \((u, v)\) would be added with multiplicity dm and weight x.

virtual_remove_edge(u, v, dm=1, entropy_args={})#

Return the difference in description length if edge \((u, v)\) with multiplicity dm would be removed.

virtual_update_edge(u, v, nx, entropy_args={})#

Return the difference in description length if edge \((u, v)\) would take a new weight nx.

virtual_update_node(v, nt, entropy_args={})#

Return the difference in description length if node v would take a new value nt.

xvals_sweep(beta=inf, niter=100, maxiter=0, tol=1e-07, min_size=1, entropy_args={})#

Perform sweeps of a greedy update on the edge weight category values, based on bisection search.

Parameters:
betafloat (optional, default: np.inf)

Inverse temperature parameter.

niterint (optional, default: 100)

Number of categories to update.

maxiterint (optional, default: 0)

Maximum number of iterations for bisection search (0 means unlimited).

tolfloat (optional, default: 1e-7)

Tolerance for bisection search.

min_sizeint (optional, default: 1)

Minimum size of edge categories that will be updated.

entropy_argsdict (optional, default: {})

Entropy arguments, with the same meaning and defaults as in entropy().

verboseboolean (optional, default: False)

If verbose == True, detailed information will be displayed.

Returns:
dSfloat

Entropy difference after the sweeps.

nmovesint

Number of variables moved.