Adjacency Matrices¶

Adjacency holds square matrices over a set of nodes: similarity and distance matrices, connectivity, directed graphs. Symmetric matrices are stored as their upper triangle and rebuilt on demand, and one object can hold a stack of matrices — one per region, subject, or timepoint. Most of its methods mirror the ones on BrainData.

Create one¶

An Adjacency accepts a numpy array, a dataframe, a CSV path, or a list of any of those. You also declare the matrix type: 'similarity' (symmetric, typically ones on the diagonal), 'distance' (symmetric, zeros on the diagonal), or 'directed' (not symmetric, stored in full). Labels are optional.

Here is fake data with three blocks of four nodes each, at signal strengths 1, 2 and 3. The noise is symmetrized before it is added, because a matrix declared symmetric has to be symmetric.

import numpy as np
from scipy.linalg import block_diag

from nltools.data import Adjacency

rng = np.random.default_rng(0)

def symmetric_noise(scale):
    noise = rng.standard_normal((12, 12)) * scale
    return (noise + noise.T) / 2

m1 = block_diag(np.ones((4, 4)), np.zeros((4, 4)), np.zeros((4, 4)))
m2 = block_diag(np.zeros((4, 4)), np.ones((4, 4)), np.zeros((4, 4)))
m3 = block_diag(np.zeros((4, 4)), np.zeros((4, 4)), np.ones((4, 4)))

blocks = Adjacency(
    (m1 * 1 + m2 * 2 + m3 * 3) + symmetric_noise(0.1),
    matrix_type="similarity",
    labels=["C1"] * 4 + ["C2"] * 4 + ["C3"] * 4,
)
blocks

nltools.data.adjacency.Adjacency(shape=(12, 12), Y=(0, 0), is_symmetric=True, matrix_type=similarity)

squareform() rebuilds the full matrix from the stored triangle:

blocks.squareform().shape

(12, 12)

plot() draws it as a heatmap, where the three blocks are visible as brighter squares on the diagonal:

blocks.plot()

cluster_summary averages edges by a grouping variable — here the labels. scope='within' summarizes edges inside each group; 'between' summarizes edges that cross groups. The three block strengths come back out.

blocks.cluster_summary(clusters=blocks.labels, scope="within", summary="mean")

{'C2': 2.0493255627514344, 'C1': 0.9698858263151181, 'C3': 3.015515137049004}

Regression¶

regress covers two different questions, and which one you get depends on what you pass as the design.

Decomposing one matrix¶

Pass an Adjacency and the matrix is decomposed into a weighted sum of the predictor matrices; edges are the observations and matrices the predictors. Using the three blocks as predictors recovers the weights the data was built with. Predictor and response must agree on node ordering, so the design carries the same labels.

design = Adjacency([m1, m2, m3], matrix_type="similarity", labels=blocks.labels)
block_fit = blocks.regress(design)

print(f"beta:  {block_fit['beta'].round(3)}")
print(f"t:     {block_fit['t'].round(1)}")
print(f"df:    {block_fit['df']}")

beta:  [0.97  2.049 3.016]
t:     [ 33.   69.6 102.5]
df:    63

Regression at every edge¶

Pass a DesignMatrix instead and each edge gets its own regression across a stack of matrices — the analogue of a mass-univariate imaging analysis, with the same multiple-comparisons problem.

The data here is 15 matrices: five with the first block on, five with it off, five on again.

from nltools.data import DesignMatrix

stack = Adjacency(
    [m1 + symmetric_noise(0.5) for _ in range(5)]
    + [symmetric_noise(0.5) for _ in range(5)]
    + [m1 + symmetric_noise(0.5) for _ in range(5)],
    matrix_type="similarity",
)

on_off_on = DesignMatrix(
    np.array([1] * 5 + [0] * 5 + [1] * 5).reshape(-1, 1),
    columns=["on"],
    sampling_freq=1.0,
)
on_off_on.plot(title="Model")

The result is one map per predictor. Plotting the t map above a cutoff shows the edges that follow the on-off-on pattern, which is the block we planted:

edge_fit = stack.regress(on_off_on)
edge_fit["t"].plot(vmin=2)

Similarity and distance¶

similarity compares two matrices and tests the result by permutation. It returns the correlation and a p-value; pass random_state to make the permutations reproducible. n_jobs=1 keeps the 5,000 permutations in this process: on a matrix this small, starting worker processes costs more than the parallelism saves.

blocks.similarity(
    Adjacency(m1, matrix_type="similarity", labels=blocks.labels),
    metric="spearman",
    n_jobs=1,
    random_state=0,
)

{'correlation': np.float64(0.29879205721671714), 'p': 0.007598480303939212}

distance measures every matrix in a stack against every other, returning a new Adjacency whose nodes are the original matrices. Any metric scikit-learn's pairwise_distances accepts works. The five on-matrices, the five off-matrices and the final five on-matrices group as you would hope.

matrix_distance = stack.distance(metric="correlation")
matrix_distance.plot()

distance_to_similarity converts back the other way:

matrix_distance.distance_to_similarity(metric="correlation").plot()

Multidimensional scaling¶

A distance matrix can be laid out in space. plot_mds does that in two or three dimensions, which is a quick way to see whether the on and off matrices separate.

labeled_distance = matrix_distance.copy()
labeled_distance.labels = ["On"] * 5 + ["Off"] * 5 + ["On"] * 5
labeled_distance.plot_mds(n_components=3)

Graphs¶

to_graph hands a matrix to networkx, so every graph metric and layout is available. Here the three noiseless blocks make three disconnected cliques of four nodes, and every node has degree 3.

import networkx as nx

clique_graph = Adjacency(m1 + m2 + m3, matrix_type="similarity").to_graph()
print(f"degree of each node: {dict(clique_graph.degree())}")

nx.draw_circular(clique_graph, node_color="lightsteelblue", with_labels=True)

degree of each node: {0: 3, 1: 3, 2: 3, 3: 3, 4: 3, 5: 3, 6: 3, 7: 3, 8: 3, 9: 3, 10: 3, 11: 3}