mne_connectivity.learn_graph#

mne_connectivity.learn_graph(data, *, names=None, indices='lower', alpha=0.5, beta=0.5, metric='sqeuclidean', normalize=True, max_iter=1000, tol=0.0001, step_size=None, n_jobs=1, verbose=None)[source]#

Learn a sparse graph from smooth signals.

Estimates a symmetric, non-negative, zero-diagonal weighted adjacency matrix directly from signal smoothness, rather than from pairwise similarity/coherence computed independently per edge. All edge weights are learned jointly by solving

\[\min_{w \geq 0} \; 2 w^\top z - \alpha \mathbf{1}^\top \log(Sw) + \beta \|w\|^2\]

where \(z\) is the vector of pairwise distances between nodes (rows of data) and \(S\) maps an edge-weight vector to the corresponding vector of node degrees. This is the model proposed in [1], and is the same model used on real MEG data (with \(\alpha\), \(\beta\) named \(\sigma\), \(\rho\)) in [2], where it was found to be more robust than Phase Lag Index for short recordings and to better capture individual-specific connectivity. The default alpha=0.5, beta=0.5 is one of the parameter combinations swept in the robustness analysis of [2]. See Notes for details.

Parameters:
dataarray_like, shape ([n_epochs,] n_nodes, n_times) | Epochs | generator

The data from which to learn a graph. Can be a single 2D array of shape (n_nodes, n_times), a 3D array of shape (n_epochs, n_nodes, n_times), a mne.Epochs object, or a list/generator of 2D arrays or mne.SourceEstimate / mne.VolSourceEstimate objects. For mne.Epochs, channels marked as bad are not used to learn the graph; they are kept as nodes, with their connections set to NaN.

namesarray_like | None

The names of the nodes of the dataset used to compute connectivity. If None (default), then names will be a list of integers from 0 to n_nodes. If a list of names, then it must be equal in length to n_nodes. Ignored (and overridden by data.ch_names) if data is an mne.Epochs object.

indices'lower' | 'upper' | 'all'

Which part of the symmetric graph to store. If 'lower' or 'upper', only the lower- or upper-triangular part is stored, and the rest is filled in when requesting get_data(output='dense'). If 'all', the full matrix (including the zero diagonal) is stored. The graph itself is always learned over all pairs of nodes. Default is 'lower'.

alphafloat

Controls how strongly node degrees are pulled away from zero via the log-barrier term. Must be strictly positive. Larger values yield more connected (higher-degree) graphs; every node is guaranteed at least one edge for any alpha > 0. Default 0.5.

betafloat

Controls the \(\|w\|^2\) penalty on edge weights. Must be non-negative. Larger values yield denser graphs. beta=0 yields the sparsest graph consistent with keeping every node connected, but lacks a strongly-convex term and so converges markedly more slowly (see max_iter); the default beta=0.5 gives fast, reliable convergence while still yielding a moderately sparse graph. See Notes.

metricstr | callable()

Distance metric passed to scipy.spatial.distance.pdist() to build the pairwise dissimilarity vector between nodes. Must yield only non-negative values. Default 'sqeuclidean' (squared Euclidean distance), matching [1] and [2].

normalizebool

If True (default), rescale the pairwise-distance vector by its mean before solving, so that the default alpha/beta behave consistently regardless of the input’s physical units or amplitude scale. See Notes.

max_iterint

Maximum number of primal-dual iterations. Default 1000. Sparser settings (small beta, in particular beta=0) may require more iterations to converge. See Notes.

tolfloat

Relative convergence tolerance on both the edge-weight and node-degree updates. Default 1e-4.

step_sizefloat | None

Override the closed-form step size 1 / (1 + 2 * beta + sqrt(2 * (n_nodes - 1))). Must be strictly positive. Useful for matching the step size used by another implementation (e.g. for validation) or for trading off convergence speed against stability. Default None.

n_jobsint

The number of jobs to run in parallel (default 1). Requires the joblib package. Each epoch’s graph is learned independently, so this parallelizes across epochs; for small numbers of nodes/iterations the joblib dispatch overhead may outweigh the benefit.

verbosebool | str | int | None

If not None, override default verbose level (see mne.verbose() for more info). If used, it should be passed as a keyword-argument only.

Returns:
graphinstance of Connectivity or EpochConnectivity

The learned graph(s): a symmetric, zero-diagonal, non-negative weighted adjacency, stored as the raveled values of the part of the matrix selected by indices. A single Connectivity is returned if data was a bare 2D array; otherwise an EpochConnectivity is returned with one graph per epoch.

See also

mne_connectivity.Connectivity
mne_connectivity.EpochConnectivity

Notes

This implements the model of eq. 12 in [1] (their \(\alpha\), \(\beta\)), identical to the model used in [2] (that paper’s \(\sigma\), \(\rho\)). The log-barrier term on node degrees (rather than on individual edge weights) is what distinguishes this model from a plain Gaussian-kernel graph construction: it guarantees every node has at least one edge, without preventing individual edges from being exactly zero, so sparse and fully connected graphs can be learned simultaneously.

Because \(z\) (the pairwise-distance vector) scales with the physical units/amplitude of data, and alpha/beta interact multiplicatively with that scale, normalize=True divides \(z\) by its mean before solving so that the default alpha=0.5, beta=0.5 are reasonable starting points regardless of the input’s scale. Use normalize=False for a literal, unscaled port of the optimization problem as stated above.

When beta=0, the learned graph is the sparsest one consistent with keeping every node connected, and changing alpha alone only rescales the overall magnitude of the learned weights, not the sparsity pattern (see Proposition 2 of [1]). However, beta=0 removes the only strongly-convex term in the objective, so this primal-dual solver converges much more slowly in that regime (often thousands of iterations) and the relative-change stopping criterion (tol) is a weaker guarantee of closeness to the true optimum than it is for beta > 0. For most uses, prefer sweeping beta away from 0 (as the default does) for both faster and more reliable convergence.

The optimization is solved per epoch using a forward-backward-forward primal-dual splitting with per-iteration complexity O(n_nodes**2).

References

Examples using mne_connectivity.learn_graph#

Learn a functional connectivity graph

Learn a functional connectivity graph