algorithms.graph.field¶
Module: algorithms.graph.field
¶
Inheritance diagram for nipy.algorithms.graph.field
:
This module implements the Field class, which simply a WeightedGraph (see the graph.py) module, plus an array that yields (possibly multi-dimnesional) features associated with graph vertices. This allows some kinds of computations (all those relating to mathematical morphology, diffusion etc.)
Certain functions are provided to Instantiate Fields easily, given a WeightedGraph and feature data.
Author:Bertrand Thirion, 2006–2011
Class¶
Field
¶
- class nipy.algorithms.graph.field.Field(V, edges=None, weights=None, field=None)¶
Bases:
WeightedGraph
- This is the basic field structure,
which contains the weighted graph structure plus an array of data (the ‘field’)
- field is an array of size(n, p)
where n is the number of vertices of the graph and p is the field dimension
- __init__(V, edges=None, weights=None, field=None)¶
- Parameters:
- V (int > 0) the number of vertices of the graph
- edges=None: the edge array of the graph
- weights=None: the associated weights array
- field=None: the field data itself
- adjacency()¶
returns the adjacency matrix of the graph as a sparse coo matrix
- Returns:
- adj: scipy.sparse matrix instance,
that encodes the adjacency matrix of self
- anti_symmeterize()¶
anti-symmeterize self, i.e. produces the graph whose adjacency matrix would be the antisymmetric part of its current adjacency matrix
- cc()¶
Compte the different connected components of the graph.
- Returns:
- label: array of shape(self.V), labelling of the vertices
- cliques()¶
Extraction of the graphe cliques these are defined using replicator dynamics equations
- Returns:
- cliques: array of shape (self.V), type (np.int_)
labelling of the vertices according to the clique they belong to
- closing(nbiter=1)¶
Morphological closing of the field data. self.field is changed inplace
- Parameters:
- nbiter=1the number of iterations required
- compact_neighb()¶
returns a compact representation of self
- Returns:
- idx: array of of shape(self.V + 1):
the positions where to find the neighbors of each node within neighb and weights
- neighb: array of shape(self.E), concatenated list of neighbors
- weights: array of shape(self.E), concatenated list of weights
- constrained_voronoi(seed)¶
Voronoi parcellation of the field starting from the input seed
- Parameters:
- seed: int array of shape(p), the input seeds
- Returns:
- label: The resulting labelling of the data
Notes
FIXME: deal with graphs with several ccs
- copy()¶
copy function
- custom_watershed(refdim=0, th=-inf)¶
customized watershed analysis of the field. Note that bassins are found around each maximum (and not minimum as conventionally)
- Parameters:
- refdim: int, optional
- th: float optional, threshold of the field
- Returns:
- idx: array of shape (nbassins)
indices of the vertices that are local maxima
- labelarray of shape (self.V)
labelling of the vertices according to their bassin
- cut_redundancies()¶
Returns a graph with redundant edges removed: ecah edge (ab) is present only once in the edge matrix: the correspondng weights are added.
- Returns:
- the resulting WeightedGraph
- degrees()¶
Returns the degree of the graph vertices.
- Returns:
- rdegree: (array, type=int, shape=(self.V,)), the right degrees
- ldegree: (array, type=int, shape=(self.V,)), the left degrees
- diffusion(nbiter=1)¶
diffusion of the field data in the weighted graph structure self.field is changed inplace
- Parameters:
- nbiter: int, optional the number of iterations required
Notes
The process is run for all the dimensions of the field
- dijkstra(seed=0)¶
Returns all the [graph] geodesic distances starting from seed x
- seed (int, >-1, <self.V) or array of shape(p)
edge(s) from which the distances are computed
- Returns:
- dg: array of shape (self.V),
the graph distance dg from ant vertex to the nearest seed
Notes
It is mandatory that the graph weights are non-negative
- dilation(nbiter=1, fast=True)¶
Morphological dilation of the field data, changed in place
- Parameters:
- nbiter: int, optional, the number of iterations required
Notes
When data dtype is not float64, a slow version of the code is used
- erosion(nbiter=1)¶
Morphological opening of the field
- Parameters:
- nbiter: int, optional, the number of iterations required
- floyd(seed=None)¶
Compute all the geodesic distances starting from seeds
- Parameters:
- seed= None: array of shape (nbseed), type np.int_
vertex indexes from which the distances are computed if seed==None, then every edge is a seed point
- Returns:
- dg array of shape (nbseed, self.V)
the graph distance dg from each seed to any vertex
Notes
It is mandatory that the graph weights are non-negative. The algorithm proceeds by repeating Dijkstra’s algo for each seed. Floyd’s algo is not used (O(self.V)^3 complexity…)
- from_3d_grid(xyz, k=18)¶
Sets the graph to be the topological neighbours graph of the three-dimensional coordinates set xyz, in the k-connectivity scheme
- Parameters:
- xyz: array of shape (self.V, 3) and type np.int_,
- k = 18: the number of neighbours considered. (6, 18 or 26)
- Returns:
- E(int): the number of edges of self
- geodesic_kmeans(seeds=None, label=None, maxiter=100, eps=0.0001, verbose=0)¶
Geodesic k-means algorithm i.e. obtention of clusters that are topologically connected and minimally variable concerning the information of self.field
- Parameters:
- seeds: array of shape(p), optional,
initial indices of the seeds within the field if seeds==None the labels are used as initialization
- labels: array of shape(self.V) initial labels, optional,
it is expected that labels take their values in a certain range (0..lmax) if Labels==None, this is not used if seeds==None and labels==None, an ewxception is raised
- maxiter: int, optional,
maximal number of iterations
- eps: float, optional,
increase of inertia at which convergence is declared
- Returns:
- seeds: array of shape (p), the final seeds
- labelarray of shape (self.V), the resulting field label
- J: float, inertia value
- get_E()¶
To get the number of edges in the graph
- get_V()¶
To get the number of vertices in the graph
- get_edges()¶
To get the graph’s edges
- get_field()¶
- get_local_maxima(refdim=0, th=-inf)¶
Look for the local maxima of one dimension (refdim) of self.field
- Parameters:
- refdim (int) the field dimension over which the maxima are looked after
- th = float, optional
threshold so that only values above th are considered
- Returns:
- idx: array of shape (nmax)
indices of the vertices that are local maxima
- depth: array of shape (nmax)
topological depth of the local maxima : depth[idx[i]] = q means that idx[i] is a q-order maximum
- get_vertices()¶
To get the graph’s vertices (as id)
- get_weights()¶
- highest_neighbor(refdim=0)¶
Computes the neighbor with highest field value along refdim
- Parameters:
- refdim: int, optional,
the dimension of the field under consideration
- Returns:
- hneighb: array of shape(self.V),
index of the neighbor with highest value
- is_connected()¶
States whether self is connected or not
- kruskal()¶
Creates the Minimum Spanning Tree of self using Kruskal’s algo. efficient is self is sparse
- Returns:
- K, WeightedGraph instance: the resulting MST
Notes
If self contains several connected components, will have the same number k of connected components
- left_incidence()¶
Return left incidence matrix
- Returns:
- left_incid: list
the left incidence matrix of self as a list of lists: i.e. the list[[e.0.0, .., e.0.i(0)], .., [e.V.0, E.V.i(V)]] where e.i.j is the set of edge indexes so that e.i.j[0] = i
- list_of_neighbors()¶
returns the set of neighbors of self as a list of arrays
- local_maxima(refdim=0, th=-inf)¶
Returns all the local maxima of a field
- Parameters:
- refdim (int) field dimension over which the maxima are looked after
- th: float, optional
threshold so that only values above th are considered
- Returns:
- depth: array of shape (nmax)
a labelling of the vertices such that depth[v] = 0 if v is not a local maximum depth[v] = 1 if v is a first order maximum … depth[v] = q if v is a q-order maximum
- main_cc()¶
Returns the indexes of the vertices within the main cc
- Returns:
- idx: array of shape (sizeof main cc)
- normalize(c=0)¶
Normalize the graph according to the index c Normalization means that the sum of the edges values that go into or out each vertex must sum to 1
- Parameters:
- c=0 in {0, 1, 2}, optional: index that designates the way
according to which D is normalized c == 0 => for each vertex a, sum{edge[e, 0]=a} D[e]=1 c == 1 => for each vertex b, sum{edge[e, 1]=b} D[e]=1 c == 2 => symmetric (‘l2’) normalization
Notes
Note that when sum_{edge[e, .] == a } D[e] = 0, nothing is performed
- opening(nbiter=1)¶
Morphological opening of the field data. self.field is changed inplace
- Parameters:
- nbiter: int, optional, the number of iterations required
- remove_edges(valid)¶
Removes all the edges for which valid==0
- Parameters:
- valid(self.E,) array
- remove_trivial_edges()¶
Removes trivial edges, i.e. edges that are (vv)-like self.weights and self.E are corrected accordingly
- Returns:
- self.E (int): The number of edges
- right_incidence()¶
Return right incidence matrix
- Returns:
- right_incid: list
the right incidence matrix of self as a list of lists: i.e. the list[[e.0.0, .., e.0.i(0)], .., [e.V.0, E.V.i(V)]] where e.i.j is the set of edge indexes so that e.i.j[1] = i
- set_edges(edges)¶
Sets the graph’s edges
Preconditions:
edges has a correct size
edges take values in [1..V]
- set_euclidian(X)¶
Compute the weights of the graph as the distances between the corresponding rows of X, which represents an embedding of self
- Parameters:
- X array of shape (self.V, edim),
the coordinate matrix of the embedding
- set_field(field)¶
- set_gaussian(X, sigma=0)¶
Compute the weights of the graph as a gaussian function of the distance between the corresponding rows of X, which represents an embedding of self
- Parameters:
- X array of shape (self.V, dim)
the coordinate matrix of the embedding
- sigma=0, float: the parameter of the gaussian function
Notes
When sigma == 0, the following value is used:
sigma = sqrt(mean(||X[self.edges[:, 0], :]-X[self.edges[:, 1], :]||^2))
- set_weights(weights)¶
Set edge weights
- Parameters:
- weights: array
array shape(self.V): edges weights
- show(X=None, ax=None)¶
Plots the current graph in 2D
- Parameters:
- XNone or array of shape (self.V, 2)
a set of coordinates that can be used to embed the vertices in 2D. If X.shape[1]>2, a svd reduces X for display. By default, the graph is presented on a circle
- ax: None or int, optional
ax handle
- Returns:
- ax: axis handle
Notes
This should be used only for small graphs.
- subfield(valid)¶
Returns a subfield of self, with only vertices such that valid > 0
- Parameters:
- valid: array of shape (self.V),
nonzero for vertices to be retained
- Returns:
- F: Field instance,
the desired subfield of self
Notes
The vertices are renumbered as [1..p] where p = sum(valid>0) when sum(valid) == 0 then None is returned
- subgraph(valid)¶
Creates a subgraph with the vertices for which valid>0 and with the corresponding set of edges
- Parameters:
- valid, array of shape (self.V): nonzero for vertices to be retained
- Returns:
- G, WeightedGraph instance, the desired subgraph of self
Notes
The vertices are renumbered as [1..p] where p = sum(valid>0) when sum(valid==0) then None is returned
- symmeterize()¶
Symmeterize self, modify edges and weights so that self.adjacency becomes the symmetric part of the current self.adjacency.
- threshold_bifurcations(refdim=0, th=-inf)¶
Analysis of the level sets of the field: Bifurcations are defined as changes in the topology in the level sets when the level (threshold) is varied This can been thought of as a kind of Morse analysis
- Parameters:
- th: float, optional,
threshold so that only values above th are considered
- Returns:
- idx: array of shape (nlsets)
indices of the vertices that are local maxima
- height: array of shape (nlsets)
the depth of the local maxima depth[idx[i]] = q means that idx[i] is a q-order maximum Note that this is also the diameter of the basins associated with local maxima
- parents: array of shape (nlsets)
the label of the maximum which dominates each local maximum i.e. it describes the hierarchy of the local maxima
- label: array of shape (self.V)
a labelling of thevertices according to their bassin
- to_coo_matrix()¶
Return adjacency matrix as coo sparse
- Returns:
- sp: scipy.sparse matrix instance
that encodes the adjacency matrix of self
- voronoi_diagram(seeds, samples)¶
Defines the graph as the Voronoi diagram (VD) that links the seeds. The VD is defined using the sample points.
- Parameters:
- seeds: array of shape (self.V, dim)
- samples: array of shape (nsamples, dim)
Notes
By default, the weights are a Gaussian function of the distance The implementation is not optimal
- voronoi_labelling(seed)¶
Performs a voronoi labelling of the graph
- Parameters:
- seed: array of shape (nseeds), type (np.int_),
vertices from which the cells are built
- Returns:
- labels: array of shape (self.V) the labelling of the vertices
- ward(nbcluster)¶
Ward’s clustering of self
- Parameters:
- nbcluster: int,
the number of desired clusters
- Returns:
- label: array of shape (self.V)
the resulting field label
- J (float): the resulting inertia
Functions¶
- nipy.algorithms.graph.field.field_from_coo_matrix_and_data(x, data)¶
Instantiates a weighted graph from a (sparse) coo_matrix
- Parameters:
- x: (V, V) scipy.sparse.coo_matrix instance,
the input matrix
- data: array of shape (V, dim),
the field data
- Returns:
- ifield: resulting Field instance
- nipy.algorithms.graph.field.field_from_graph_and_data(g, data)¶
Instantiate a Fieldfrom a WeightedGraph plus some feature data Parameters ———- x: (V, V) scipy.sparse.coo_matrix instance,
the input matrix
- data: array of shape (V, dim),
the field data
- Returns:
- ifield: resulting field instance