summaryrefslogtreecommitdiff
path: root/src/python/doc/rips_complex_user.rst
blob: 6048cc4e7d7d799480a19cd8c81449e6a5dc04cf (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
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
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
:orphan:

.. To get rid of WARNING: document isn't included in any toctree

Rips complex user manual
=========================
Definition
----------

====================================================================  ================================  ======================
:Authors: Clément Maria, Pawel Dlotko, Vincent Rouvreau, Marc Glisse  :Since: GUDHI 2.0.0               :License: GPL v3
====================================================================  ================================  ======================

+-------------------------------------------+----------------------------------------------------------------------+
| :doc:`rips_complex_user`                  | :doc:`rips_complex_ref`                                              |
+-------------------------------------------+----------------------------------------------------------------------+

The `Rips complex <https://en.wikipedia.org/wiki/Vietoris%E2%80%93Rips_complex>`_ is a simplicial complex that
generalizes proximity (:math:`\varepsilon`-ball) graphs to higher dimensions. The vertices correspond to the input
points, and a simplex is present if and only if its diameter is smaller than some parameter α.  Considering all
parameters α defines a filtered simplicial complex, where the filtration value of a simplex is its diameter.
The filtration can be restricted to values α smaller than some threshold, to reduce its size.  Beware that some
people define the Rips complex using a bound of 2α instead of α, particularly when comparing it to an ambient
Čech complex.  They end up with the same combinatorial object, but filtration values which are half of ours.

The input discrete metric space can be provided as a point cloud plus a distance function, or as a distance matrix.

When creating a simplicial complex from the graph, :doc:`RipsComplex <rips_complex_ref>` first builds the graph and
inserts it into the data structure. It then expands the simplicial complex (adds the simplices corresponding to cliques)
when required. The expansion can be stopped at dimension `max_dimension`, by default 1.

A vertex name corresponds to the index of the point in the given range (aka. the point cloud).

.. figure::
    ../../doc/Rips_complex/rips_complex_representation.png
    :align: center

    Rips-complex one skeleton graph representation

On this example, as edges (4,5), (4,6) and (5,6) are in the complex, simplex (4,5,6) is added with the filtration value
set with :math:`max(filtration(4,5), filtration(4,6), filtration(5,6))`. And so on for simplex (0,1,2,3).

If the :doc:`RipsComplex <rips_complex_ref>` interfaces are not detailed enough for your need, please refer to
rips_persistence_step_by_step.cpp C++ example, where the graph construction over the Simplex_tree is more detailed.

A Rips complex can easily become huge, even if we limit the length of the edges
and the dimension of the simplices. One easy trick, before building a Rips
complex on a point cloud, is to call :func:`~gudhi.sparsify_point_set` which removes points
that are too close to each other. This does not change its persistence diagram
by more than the length used to define "too close".

A more general technique is to use a sparse approximation of the Rips
introduced by Don Sheehy :cite:`sheehy13linear`. We are using the version
described in :cite:`buchet16efficient` (except that we multiply all filtration
values by 2, to match the usual Rips complex). :cite:`cavanna15geometric` proves
a :math:`\frac{1}{1-\varepsilon}`-interleaving, although in practice the
error is usually smaller.  A more intuitive presentation of the idea is
available in :cite:`cavanna15geometric`, and in a video
:cite:`cavanna15visualizing`. Passing an extra argument `sparse=0.3` at the
construction of a :class:`~gudhi.RipsComplex` object asks it to build a sparse Rips with
parameter :math:`\varepsilon=0.3`, while the default `sparse=None` builds the
regular Rips complex.


Point cloud
-----------

Example from a point cloud
^^^^^^^^^^^^^^^^^^^^^^^^^^

This example builds the neighborhood graph from the given points, up to max_edge_length.
Then it creates a :doc:`SimplexTree <simplex_tree_ref>` with it.

Finally, it is asked to display information about the simplicial complex.

.. testcode::

    import gudhi
    rips_complex = gudhi.RipsComplex(points=[[1, 1], [7, 0], [4, 6], [9, 6], [0, 14], [2, 19], [9, 17]],
                                     max_edge_length=12.0)

    simplex_tree = rips_complex.create_simplex_tree(max_dimension=1)
    result_str = 'Rips complex is of dimension ' + repr(simplex_tree.dimension()) + ' - ' + \
        repr(simplex_tree.num_simplices()) + ' simplices - ' + \
        repr(simplex_tree.num_vertices()) + ' vertices.'
    print(result_str)
    fmt = '%s -> %.2f'
    for filtered_value in simplex_tree.get_filtration():
        print(fmt % tuple(filtered_value))

When launching (Rips maximal distance between 2 points is 12.0, is expanded
until dimension 1 - one skeleton graph in other words), the output is:

.. testoutput::

    Rips complex is of dimension 1 - 18 simplices - 7 vertices.
    [0] -> 0.00
    [1] -> 0.00
    [2] -> 0.00
    [3] -> 0.00
    [4] -> 0.00
    [5] -> 0.00
    [6] -> 0.00
    [2, 3] -> 5.00
    [4, 5] -> 5.39
    [0, 2] -> 5.83
    [0, 1] -> 6.08
    [1, 3] -> 6.32
    [1, 2] -> 6.71
    [5, 6] -> 7.28
    [2, 4] -> 8.94
    [0, 3] -> 9.43
    [4, 6] -> 9.49
    [3, 6] -> 11.00

Notice that if we use

.. code-block:: python

    rips_complex = gudhi.RipsComplex(points=[[1, 1], [7, 0], [4, 6], [9, 6], [0, 14], [2, 19], [9, 17]],
                                     max_edge_length=12.0, sparse=2)

asking for a very sparse version (theory only gives some guarantee on the meaning of the output if `sparse<1`),
2 to 5 edges disappear, depending on the random vertex used to start the sparsification.

Example from OFF file
^^^^^^^^^^^^^^^^^^^^^

This example builds the :doc:`RipsComplex <rips_complex_ref>` from the given
points in an OFF file, and max_edge_length value.
Then it creates a :doc:`SimplexTree <simplex_tree_ref>` with it.

Finally, it is asked to display information about the Rips complex.


.. testcode::

    import gudhi
    off_file = gudhi.__root_source_dir__ + '/data/points/alphacomplexdoc.off'
    point_cloud = gudhi.read_points_from_off_file(off_file = off_file)
    rips_complex = gudhi.RipsComplex(points=point_cloud, max_edge_length=12.0)
    simplex_tree = rips_complex.create_simplex_tree(max_dimension=1)
    result_str = 'Rips complex is of dimension ' + repr(simplex_tree.dimension()) + ' - ' + \
        repr(simplex_tree.num_simplices()) + ' simplices - ' + \
        repr(simplex_tree.num_vertices()) + ' vertices.'
    print(result_str)
    fmt = '%s -> %.2f'
    for filtered_value in simplex_tree.get_filtration():
        print(fmt % tuple(filtered_value))

the program output is:

.. testoutput::

    Rips complex is of dimension 1 - 18 simplices - 7 vertices.
    [0] -> 0.00
    [1] -> 0.00
    [2] -> 0.00
    [3] -> 0.00
    [4] -> 0.00
    [5] -> 0.00
    [6] -> 0.00
    [2, 3] -> 5.00
    [4, 5] -> 5.39
    [0, 2] -> 5.83
    [0, 1] -> 6.08
    [1, 3] -> 6.32
    [1, 2] -> 6.71
    [5, 6] -> 7.28
    [2, 4] -> 8.94
    [0, 3] -> 9.43
    [4, 6] -> 9.49
    [3, 6] -> 11.00

Distance matrix
---------------

Example from a distance matrix
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

This example builds the one skeleton graph from the given distance matrix, and max_edge_length value.
Then it creates a :doc:`SimplexTree <simplex_tree_ref>` with it.

Finally, it is asked to display information about the simplicial complex.

.. testcode::

    import gudhi
    rips_complex = gudhi.RipsComplex(distance_matrix=[[],
                                                      [6.0827625303],
                                                      [5.8309518948, 6.7082039325],
                                                      [9.4339811321, 6.3245553203, 5],
                                                      [13.0384048104, 15.6524758425, 8.94427191, 12.0415945788],
                                                      [18.0277563773, 19.6468827044, 13.152946438, 14.7648230602, 5.3851648071],
                                                      [17.88854382, 17.1172427686, 12.0830459736, 11, 9.4868329805, 7.2801098893]],
                                     max_edge_length=12.0)

    simplex_tree = rips_complex.create_simplex_tree(max_dimension=1)
    result_str = 'Rips complex is of dimension ' + repr(simplex_tree.dimension()) + ' - ' + \
        repr(simplex_tree.num_simplices()) + ' simplices - ' + \
        repr(simplex_tree.num_vertices()) + ' vertices.'
    print(result_str)
    fmt = '%s -> %.2f'
    for filtered_value in simplex_tree.get_filtration():
        print(fmt % tuple(filtered_value))

When launching (Rips maximal distance between 2 points is 12.0, is expanded
until dimension 1 - one skeleton graph in other words), the output is:

.. testoutput::

    Rips complex is of dimension 1 - 18 simplices - 7 vertices.
    [0] -> 0.00
    [1] -> 0.00
    [2] -> 0.00
    [3] -> 0.00
    [4] -> 0.00
    [5] -> 0.00
    [6] -> 0.00
    [2, 3] -> 5.00
    [4, 5] -> 5.39
    [0, 2] -> 5.83
    [0, 1] -> 6.08
    [1, 3] -> 6.32
    [1, 2] -> 6.71
    [5, 6] -> 7.28
    [2, 4] -> 8.94
    [0, 3] -> 9.43
    [4, 6] -> 9.49
    [3, 6] -> 11.00

Example from csv file
^^^^^^^^^^^^^^^^^^^^^

This example builds the :doc:`RipsComplex <rips_complex_ref>` from the given
distance matrix in a csv file, and max_edge_length value.
Then it creates a :doc:`SimplexTree <simplex_tree_ref>` with it.

Finally, it is asked to display information about the Rips complex.


.. testcode::

    import gudhi
    distance_matrix = gudhi.read_lower_triangular_matrix_from_csv_file(csv_file=gudhi.__root_source_dir__ + \
        '/data/distance_matrix/full_square_distance_matrix.csv')
    rips_complex = gudhi.RipsComplex(distance_matrix=distance_matrix, max_edge_length=12.0)
    simplex_tree = rips_complex.create_simplex_tree(max_dimension=1)
    result_str = 'Rips complex is of dimension ' + repr(simplex_tree.dimension()) + ' - ' + \
        repr(simplex_tree.num_simplices()) + ' simplices - ' + \
        repr(simplex_tree.num_vertices()) + ' vertices.'
    print(result_str)
    fmt = '%s -> %.2f'
    for filtered_value in simplex_tree.get_filtration():
        print(fmt % tuple(filtered_value))

the program output is:

.. testoutput::

    Rips complex is of dimension 1 - 18 simplices - 7 vertices.
    [0] -> 0.00
    [1] -> 0.00
    [2] -> 0.00
    [3] -> 0.00
    [4] -> 0.00
    [5] -> 0.00
    [6] -> 0.00
    [2, 3] -> 5.00
    [4, 5] -> 5.39
    [0, 2] -> 5.83
    [0, 1] -> 6.08
    [1, 3] -> 6.32
    [1, 2] -> 6.71
    [5, 6] -> 7.28
    [2, 4] -> 8.94
    [0, 3] -> 9.43
    [4, 6] -> 9.49
    [3, 6] -> 11.00

Correlation matrix
------------------

Example from  a correlation matrix
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Analogously to the case of distance matrix, Rips complexes can be also constructed based on correlation matrix.
Given a correlation matrix M, comportment-wise 1-M is a distance matrix.
This example builds the one skeleton graph from the given corelation matrix and threshold value.
Then it creates a :doc:`SimplexTree <simplex_tree_ref>` with it.

Finally, it is asked to display information about the simplicial complex.

.. testcode::

    import gudhi
    import numpy as np

    # User defined correlation matrix is:
    # |1     0.06    0.23    0.01    0.89|
    # |0.06  1       0.74    0.01    0.61|
    # |0.23  0.74    1       0.72    0.03|
    # |0.01  0.01    0.72    1       0.7 |
    # |0.89  0.61    0.03    0.7     1   |
    correlation_matrix=np.array([[1., 0.06, 0.23, 0.01, 0.89],
                                [0.06, 1., 0.74, 0.01, 0.61],
                                [0.23, 0.74, 1., 0.72, 0.03],
                                [0.01, 0.01, 0.72, 1., 0.7],
                                [0.89, 0.61, 0.03, 0.7, 1.]], float)

    distance_matrix = 1 - correlation_matrix
    rips_complex = gudhi.RipsComplex(distance_matrix=distance_matrix, max_edge_length=1.0)

    simplex_tree = rips_complex.create_simplex_tree(max_dimension=1)
    result_str = 'Rips complex is of dimension ' + repr(simplex_tree.dimension()) + ' - ' + \
        repr(simplex_tree.num_simplices()) + ' simplices - ' + \
        repr(simplex_tree.num_vertices()) + ' vertices.'
    print(result_str)
    fmt = '%s -> %.2f'
    for filtered_value in simplex_tree.get_filtration():
        print(fmt % tuple(filtered_value))

When launching (Rips maximal distance between 2 points is 12.0, is expanded
until dimension 1 - one skeleton graph in other words), the output is:

.. testoutput::

    Rips complex is of dimension 1 - 15 simplices - 5 vertices.
    [0] -> 0.00
    [1] -> 0.00
    [2] -> 0.00
    [3] -> 0.00
    [4] -> 0.00
    [0, 4] -> 0.11
    [1, 2] -> 0.26
    [2, 3] -> 0.28
    [3, 4] -> 0.30
    [1, 4] -> 0.39
    [0, 2] -> 0.77
    [0, 1] -> 0.94
    [2, 4] -> 0.97
    [0, 3] -> 0.99
    [1, 3] -> 0.99

.. note::
    If you compute the persistence diagram and convert distances back to correlation values,
    points in the persistence diagram will be under the diagonal, and
    bottleneck distance and persistence graphical tool will not work properly,
    this is a known issue.

Weighted Rips Complex
---------------------

`WeightedRipsComplex <rips_complex_ref.html#weighted-rips-complex-reference-manual>`_ builds a simplicial complex from a distance matrix and weights on vertices.


Example from a distance matrix and weights
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

The following example computes the weighted Rips filtration associated with a distance matrix and weights on vertices.

.. testcode::

    from gudhi.weighted_rips_complex import WeightedRipsComplex
    dist = [[], [1]]
    weights = [1, 100]
    w_rips = WeightedRipsComplex(distance_matrix=dist, weights=weights)
    st = w_rips.create_simplex_tree(max_dimension=2)
    print(list(st.get_filtration()))

The output is:

.. testoutput::

    [([0], 2.0), ([1], 200.0), ([0, 1], 200.0)]

Example from a point cloud combined with DistanceToMeasure
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

Combining with DistanceToMeasure, one can compute the DTM-filtration of a point set, as in `this notebook <https://github.com/GUDHI/TDA-tutorial/blob/master/Tuto-GUDHI-DTM-filtrations.ipynb>`_. 
Remark that `DTMRipsComplex <rips_complex_user.html#dtm-rips-complex>`_ class provides exactly this function.

.. testcode::

    import numpy as np
    from scipy.spatial.distance import cdist
    from gudhi.point_cloud.dtm import DistanceToMeasure
    from gudhi.weighted_rips_complex import WeightedRipsComplex
    pts = np.array([[2.0, 2.0], [0.0, 1.0], [3.0, 4.0]])
    dist = cdist(pts,pts)
    dtm = DistanceToMeasure(2, q=2, metric="precomputed")
    r = dtm.fit_transform(dist)
    w_rips = WeightedRipsComplex(distance_matrix=dist, weights=r)
    st = w_rips.create_simplex_tree(max_dimension=2)
    print(st.persistence())

The output is:

.. testoutput::

    [(0, (3.1622776601683795, inf)), (0, (3.1622776601683795, 5.39834563766817)), (0, (3.1622776601683795, 5.39834563766817))]

DTM Rips Complex
----------------

:class:`~gudhi.dtm_rips_complex.DTMRipsComplex` builds a simplicial complex from a point set or a full distance matrix (in the form of ndarray), as described in the above example.
This class constructs a weighted Rips complex giving larger weights to outliers, which reduces their impact on the persistence diagram. See `this notebook <https://github.com/GUDHI/TDA-tutorial/blob/master/Tuto-GUDHI-DTM-filtrations.ipynb>`_ for some experiments.

.. testcode::

    import numpy as np
    from gudhi.dtm_rips_complex import DTMRipsComplex
    pts = np.array([[2.0, 2.0], [0.0, 1.0], [3.0, 4.0]])
    dtm_rips = DTMRipsComplex(points=pts, k=2)
    st = dtm_rips.create_simplex_tree(max_dimension=2)
    print(st.persistence())

The output is:

.. testoutput::

    [(0, (3.1622776601683795, inf)), (0, (3.1622776601683795, 5.39834563766817)), (0, (3.1622776601683795, 5.39834563766817))]