__init__.py 9.42 KB
Newer Older
Tiago Peixoto's avatar
Tiago Peixoto committed
1
#! /usr/bin/env python
2
# -*- coding: utf-8 -*-
Tiago Peixoto's avatar
Tiago Peixoto committed
3
#
4
5
# graph_tool -- a general graph manipulation python module
#
Tiago Peixoto's avatar
Tiago Peixoto committed
6
# Copyright (C) 2006-2013 Tiago de Paula Peixoto <tiago@skewed.de>
Tiago Peixoto's avatar
Tiago Peixoto committed
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.

"""
``graph_tool.spectral`` - Spectral properties
---------------------------------------------
24
25
26
27
28
29
30
31
32
33
34
35
36

Summary
+++++++

.. autosummary::
   :nosignatures:

   adjacency
   laplacian
   incidence

Contents
++++++++
Tiago Peixoto's avatar
Tiago Peixoto committed
37
38
"""

39
40
from __future__ import division, absolute_import, print_function

41
from .. import _degree, _prop, Graph, _limit_args
42
43
from .. stats import label_self_loops
import numpy
Tiago Peixoto's avatar
Tiago Peixoto committed
44
45
import scipy.sparse

46
47
from .. dl_import import dl_import
dl_import("from . import libgraph_tool_spectral")
Tiago Peixoto's avatar
Tiago Peixoto committed
48
49
50

__all__ = ["adjacency", "laplacian", "incidence"]

51

52
def adjacency(g, weight=None, index=None):
53
54
55
56
    r"""Return the adjacency matrix of the graph.

    Parameters
    ----------
57
    g : :class:`~graph_tool.Graph`
58
        Graph to be used.
59
    weight : :class:`~graph_tool.PropertyMap` (optional, default: True)
60
        Edge property map with the edge weights.
61
62
63
    index : :class:`~graph_tool.PropertyMap` (optional, default: None)
        Vertex property map specifying the row/column indexes. If not provided, the
        internal vertex index is used.
64
65
66

    Returns
    -------
67
68
    a : :mod:`~scipy.sparse.csr_matrix`
        The (sparse) adjacency matrix.
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84

    Notes
    -----
    The adjacency matrix is defined as

    .. math::

        a_{i,j} =
        \begin{cases}
            1 & \text{if } v_i \text{ is adjacent to } v_j, \\
            0 & \text{otherwise}
        \end{cases}

    In the case of weighted edges, the value 1 is replaced the weight of the
    respective edge.

85
86
87
    In the case of networks with parallel edges, the entries in the matrix
    become simply the edge multiplicities.

88
89
    Examples
    --------
90
91
92
93
94
    .. testsetup::

      gt.seed_rng(42)

    >>> g = gt.random_graph(100, lambda: (10, 10))
95
    >>> m = gt.adjacency(g)
96
    >>> print(m.todense())
Tiago Peixoto's avatar
Tiago Peixoto committed
97
    [[ 0.  0.  0. ...,  0.  1.  0.]
98
     [ 0.  0.  0. ...,  0.  0.  0.]
Tiago Peixoto's avatar
Tiago Peixoto committed
99
     [ 0.  0.  0. ...,  0.  0.  1.]
Tiago Peixoto's avatar
Tiago Peixoto committed
100
     ..., 
101
     [ 0.  0.  0. ...,  0.  0.  0.]
Tiago Peixoto's avatar
Tiago Peixoto committed
102
103
     [ 0.  0.  0. ...,  0.  0.  0.]
     [ 0.  0.  1. ...,  0.  0.  0.]]
104
105
106

    References
    ----------
107
    .. [wikipedia-adjacency] http://en.wikipedia.org/wiki/Adjacency_matrix
108
109
    """

110
111
112
113
114
115
    if index is None:
        if g.get_vertex_filter()[0] != None:
            index = g.new_vertex_property("int64_t")
            index.fa = numpy.arange(g.num_vertices())
        else:
            index = g.vertex_index
Tiago Peixoto's avatar
Tiago Peixoto committed
116

117
118
119
120
    E = g.num_edges() if g.is_directed() else 2 * g.num_edges()
    data = numpy.zeros(E, dtype="double")
    i = numpy.zeros(E, dtype="int32")
    j = numpy.zeros(E, dtype="int32")
Tiago Peixoto's avatar
Tiago Peixoto committed
121

122
123
124
125
126
    libgraph_tool_spectral.adjacency(g._Graph__graph, _prop("v", g, index),
                                     _prop("e", g, weight), data, i, j)
    m = scipy.sparse.coo_matrix((data, (i,j)))
    m = m.tocsr()
    return m
Tiago Peixoto's avatar
Tiago Peixoto committed
127

Tiago Peixoto's avatar
Tiago Peixoto committed
128
129

@_limit_args({"deg": ["total", "in", "out"]})
130
def laplacian(g, deg="total", normalized=True, weight=None, index=None):
131
132
133
134
    r"""Return the Laplacian matrix of the graph.

    Parameters
    ----------
135
    g : :class:`~graph_tool.Graph`
136
137
138
139
140
        Graph to be used.
    deg : str (optional, default: "total")
        Degree to be used, in case of a directed graph.
    normalized : bool (optional, default: True)
        Whether to compute the normalized Laplacian.
141
    weight : :class:`~graph_tool.PropertyMap` (optional, default: True)
142
        Edge property map with the edge weights.
143
144
145
    index : :class:`~graph_tool.PropertyMap` (optional, default: None)
        Vertex property map specifying the row/column indexes. If not provided, the
        internal vertex index is used.
146
147
148

    Returns
    -------
149
150
    l : :mod:`~scipy.sparse.csr_matrix`
        The (sparse) Laplacian matrix.
151
152
153

    Notes
    -----
Tiago Peixoto's avatar
Tiago Peixoto committed
154
    The weighted Laplacian matrix is defined as
155
156
157

    .. math::

Tiago Peixoto's avatar
Tiago Peixoto committed
158
        \ell_{ij} =
159
160
        \begin{cases}
        \Gamma(v_i) & \text{if } i = j \\
Tiago Peixoto's avatar
Tiago Peixoto committed
161
        -w_{ij}     & \text{if } i \neq j \text{ and } v_i \text{ is adjacent to } v_j \\
162
163
164
        0           & \text{otherwise}.
        \end{cases}

165
166
    Where :math:`\Gamma(v_i)=\sum_j A_{ij}w_{ij}` is sum of the weights of the
    edges incident on vertex :math:`v_i`. The normalized version is
167
168
169

    .. math::

Tiago Peixoto's avatar
Tiago Peixoto committed
170
        \ell_{ij} =
171
172
        \begin{cases}
        1         & \text{ if } i = j \text{ and } \Gamma(v_i) \neq 0 \\
Tiago Peixoto's avatar
Tiago Peixoto committed
173
        -\frac{w_{ij}}{\sqrt{\Gamma(v_i)\Gamma(v_j)}} & \text{ if } i \neq j \text{ and } v_i \text{ is adjacent to } v_j \\
174
175
176
        0         & \text{otherwise}.
        \end{cases}

Tiago Peixoto's avatar
Tiago Peixoto committed
177
178
179
180
181
    In the case of unweighted edges, it is assumed :math:`w_{ij} = 1`.

    For directed graphs, it is assumed :math:`\Gamma(v_i)=\sum_j A_{ij}w_{ij} +
    \sum_j A_{ji}w_{ji}` if ``deg=="total"``, :math:`\Gamma(v_i)=\sum_j A_{ij}w_{ij}`
    if ``deg=="out"`` or :math:`\Gamma(v_i)=\sum_j A_{ji}w_{ji}` ``deg=="in"``.
182
183
184

    Examples
    --------
185
186
187
188
    .. testsetup::

      gt.seed_rng(42)

189
190
    >>> g = gt.random_graph(100, lambda: (10,10))
    >>> m = gt.laplacian(g)
191
    >>> print(m.todense())
Tiago Peixoto's avatar
Tiago Peixoto committed
192
193
194
    [[ 1.   -0.05  0.   ...,  0.    0.    0.  ]
     [ 0.    1.    0.   ...,  0.    0.   -0.05]
     [ 0.    0.    1.   ...,  0.   -0.05  0.  ]
Tiago Peixoto's avatar
Tiago Peixoto committed
195
     ..., 
Tiago Peixoto's avatar
Tiago Peixoto committed
196
197
198
     [ 0.    0.    0.   ...,  1.    0.    0.  ]
     [-0.05  0.    0.   ...,  0.    1.    0.  ]
     [ 0.    0.    0.   ..., -0.05  0.    1.  ]]
199
200
201

    References
    ----------
202
    .. [wikipedia-laplacian] http://en.wikipedia.org/wiki/Laplacian_matrix
203
204
    """

205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
    if index is None:
        if g.get_vertex_filter()[0] != None:
            index = g.new_vertex_property("int64_t")
            index.fa = numpy.arange(g.num_vertices())
        else:
            index = g.vertex_index

    nself = label_self_loops(g, mark_only=True).a.sum()
    E = g.num_edges() - nself
    if not g.is_directed():
        E *= 2
    N = E + g.num_vertices()
    data = numpy.zeros(N, dtype="double")
    i = numpy.zeros(N, dtype="int32")
    j = numpy.zeros(N, dtype="int32")

    if normalized:
        libgraph_tool_spectral.norm_laplacian(g._Graph__graph, _prop("v", g, index),
                                              _prop("e", g, weight), deg, data, i, j)
Tiago Peixoto's avatar
Tiago Peixoto committed
224
    else:
225
226
227
228
        libgraph_tool_spectral.laplacian(g._Graph__graph, _prop("v", g, index),
                                         _prop("e", g, weight), deg, data, i, j)
    m = scipy.sparse.coo_matrix((data, (i,j)))
    m = m.tocsr()
Tiago Peixoto's avatar
Tiago Peixoto committed
229
230
    return m

Tiago Peixoto's avatar
Tiago Peixoto committed
231

232
def incidence(g, vindex=None, eindex=None):
233
234
235
236
    r"""Return the incidence matrix of the graph.

    Parameters
    ----------
237
    g : :class:`~graph_tool.Graph`
238
        Graph to be used.
239
240
241
242
243
244
    vindex : :class:`~graph_tool.PropertyMap` (optional, default: None)
        Vertex property map specifying the row indexes. If not provided, the
        internal vertex index is used.
    eindex : :class:`~graph_tool.PropertyMap` (optional, default: None)
        Edge property map specifying the column indexes. If not provided, the
        internal edge index is used.
245
246
247

    Returns
    -------
248
249
    a : :mod:`~scipy.sparse.csr_matrix`
        The (sparse) adjacency matrix.
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275

    Notes
    -----
    For undirected graphs, the incidence matrix is defined as

    .. math::

        b_{i,j} =
        \begin{cases}
            1 & \text{if vertex } v_i \text{and edge } e_j \text{ are incident}, \\
            0 & \text{otherwise}
        \end{cases}

    For directed graphs, the definition is

    .. math::

        b_{i,j} =
        \begin{cases}
            1 & \text{if edge } e_j \text{ enters vertex } v_i, \\
            -1 & \text{if edge } e_j \text{ leaves vertex } v_i, \\
            0 & \text{otherwise}
        \end{cases}

    Examples
    --------
276
277
278
279
    .. testsetup::

      gt.seed_rng(42)

280
281
    >>> g = gt.random_graph(100, lambda: (2,2))
    >>> m = gt.incidence(g)
282
    >>> print(m.todense())
Tiago Peixoto's avatar
Tiago Peixoto committed
283
284
    [[-1. -1.  0. ...,  0.  0.  0.]
     [ 0.  0. -1. ...,  0.  0.  0.]
285
     [ 0.  0.  0. ...,  0.  0.  0.]
Tiago Peixoto's avatar
Tiago Peixoto committed
286
     ..., 
287
     [ 0.  0.  0. ...,  0.  0.  0.]
Tiago Peixoto's avatar
Tiago Peixoto committed
288
289
     [ 0.  0.  0. ..., -1.  0.  0.]
     [ 0.  0.  0. ...,  0. -1. -1.]]
290
291
292

    References
    ----------
293
    .. [wikipedia-incidence] http://en.wikipedia.org/wiki/Incidence_matrix
294
295
    """

296
297
298
299
300
301
    if vindex is None:
        if g.get_edge_filter()[0] != None:
            vindex = g.new_vertex_property("int64_t")
            vindex.fa = numpy.arange(g.num_vertices())
        else:
            vindex = g.vertex_index
Tiago Peixoto's avatar
Tiago Peixoto committed
302

303
304
305
306
307
308
    if eindex is None:
        if g.get_edge_filter()[0] != None:
            eindex = g.new_edge_property("int64_t")
            eindex.fa = numpy.arange(g.num_edges())
        else:
            eindex = g.edge_index
Tiago Peixoto's avatar
Tiago Peixoto committed
309
310

    E = g.num_edges()
311
312
313
314
315
316
317
318
    data = numpy.zeros(2 * E, dtype="double")
    i = numpy.zeros(2 * E, dtype="int32")
    j = numpy.zeros(2 * E, dtype="int32")

    libgraph_tool_spectral.incidence(g._Graph__graph, _prop("v", g, vindex),
                                     _prop("e", g, eindex), data, i, j)
    m = scipy.sparse.coo_matrix((data, (i,j)))
    m = m.tocsr()
Tiago Peixoto's avatar
Tiago Peixoto committed
319
    return m