Search code examples
pythonnumpyeigenvector

Changing eigenvalue normalization in numpy


I would like to write a function that returns the connected components of an undirected graph. For example, for a graph with five nodes and two connected components (nodes 0, 1, 3, 4 being connected and node 2 being isolated) it would return two vectors

[1, 1, 0, 1, 1]
[0, 0, 1, 0, 0]

with ones if two nodes are connected and zeros otherwise.

I'm using the graph Laplacian for this, which for undirected graphs is a symmetric real matrix. I'm computing its eigenvalues and eigenvectors using Numpy.

[evalues, evectors] = numpy.linalg.eigh(laplacian)
indices = evalues < 0.00000001
evectors = numpy.transpose(evectors)
relevant_evectors = evectors[indices]

eigh returns normalized eigenvectors. For example, in the above case, the relvant_evectors are:

[[-0.49316146 -0.49316146  0.16482443 -0.49316146 -0.49316146]
 [-0.08241221 -0.08241221 -0.98632292 -0.08241221 -0.08241221]]

So as a next step I want to transform the relevant_eigenvectos to contain only zeros and ones like above. I can do this by brutally searching and replacing values. But is there a more elegant way to do this?


Solution

  • An array of vector with two kinds of elements can be converted to an array of 0-1 vector with the following:

    components = np.where(relevant_evectors < np.mean(relevant_evectors, axis=1), 0, 1)
    

    Here, the mean is taken of each row (since you transposed to get row vectors). The elements below the mean are replaced with 0, the others with 1.

    That said, I could not reproduce the kind of vectors you are getting; the eigenvectors I get are all sparse, with 0 values outside a column. Perhaps I haven't considered a complex enough graph. Here is a my version of the code for completeness.

    import numpy as np
    adjacency = np.zeros((8, 8))
    edges = [(0, 1), (0, 7), (2, 3), (4, 5), (5, 6), (4, 6)]
    for edge in edges:
        adjacency[edge] = 1
    adjacency += adjacency.T      # symmetric adjacency matrix 
    laplacian = np.diag(adjacency.sum(axis=0)) - adjacency
    [evalues, evectors] = np.linalg.eigh(laplacian)
    relevant_evectors = evectors[:, evalues < 1e-7]
    components = np.where(np.abs(relevant_evectors) < 1e-7, 0, 1)  # Works in my example...
    print(components)