# Importing libraries
import giotto as go
import giotto.time_series as ts
import giotto.diagrams as diag
import giotto.homology as hl
from giotto.diagrams import PersistenceEntropy, BettiCurve, PersistenceLandscape, HeatKernel
import numpy as np
import sklearn as sk
from sklearn.pipeline import Pipeline
from sklearn.metrics import pairwise_distances
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
The plotting.py file is required to use the following plotting functions. It can be found in the examples folder on out github https://github.com/giotto-learn/giotto-learn.
# Plotting functions
from plotting import plot_diagram, plot_landscapes
from plotting import plot_betti_surfaces, plot_betti_curves
from plotting import plot_point_cloud
We are going to consider three classical topological spaces: the circle, the 2-torus and the 2-sphere. The purpose of this tutorial is to go thgough the most famous topological spaces and compute their homology groups.
Each of the topological spaces we are going to encounter will be sampled. The resulting point clound will be the input of the persistent homology pipeline. The first step is to apply the Vietoris-Rips technique to the point cloud. Finally, the persistent homology groups will be computed.
# Representing the circle in 3d with parametric equations.
circle = np.asarray([[np.sin(t),np.cos(t),0] for t in range(400)])
plot_point_cloud(circle)
# Representing the sphere in 3d with parametric equations
sphere = np.asarray([[np.cos(s)*np.cos(t),np.cos(s)*np.sin(t),np.sin(s)] for t in range(20) for s in range(20)])
plot_point_cloud(sphere)
# Representing the torus in 3d with parametric equations
torus = np.asarray([[(2+np.cos(s))*np.cos(t),(2+np.cos(s))*np.sin(t),np.sin(s)] for t in range(20) for s in range(20)])
plot_point_cloud(torus)
# Saving the results into an array
topological_spaces=np.asarray([circle,sphere,torus])
In the next section we are using giotto to compute the persistent homology groups of the topological spaces we just constructed
# the homology ranks we choose to consider
homologyDimensions = (0, 1 ,2)
persistenceDiagram = hl.VietorisRipsPersistence(metric='euclidean', max_edge_length=10,
homology_dimensions=homologyDimensions,
n_jobs=-1)
persistenceDiagram.fit(topological_spaces)
# List of all the time-pordered persistent diagrams obtained from the list of correlation matrices
Diagrams = persistenceDiagram.transform(topological_spaces)
print(Diagrams.shape)
The topological information of the point cloud is synthesised in the persistent diagram. The horizonral axis corresponds to the moment in which an homological generator is born, while the vertical axis corresponds to the moments in which an homological generator dies. The generators of the homology groups (at given rank) are colored differently
# plotting the persistent diagram of the circle
plot_diagram(Diagrams[0])
# plotting the persistent diagram of the sphere
plot_diagram(Diagrams[1])
# plotting the persistent diagram of the torus
plot_diagram(Diagrams[2])
As you can see from the persistence diagrams, all the betti numbers were found. Some other persistent generators are also appearing, depending on how dense the sampling is and how much noise there is. For example, we see a rahter neat persistent diagram over the Torus bottle (wesee 2 persistent generators for $H^1$ and 1 persistent generator for $H^2$). Notice though that there are other persistent H1 generators, possibly due to the non-uniform sampling method we used for the torus. On the other hand, the persistent diagram for the circle is as perfect as it could be: one unique generator of $H^1$ and no other persistent generator, as expeceted.
We are going to consider different classical shapes: the real projective space and the Klein bottle. The purpose of the second part of the tutorial is to define shapes via a distance matrix. We also add noise to the distace matrix: the main reason is not to have overlapping points in the persistent diagram.
Each of the topological spaces we are going to encounter will be sampled discretely. Aftewards the Vietoris-Rips technique will be applied to the surface and the persistent homology groups will be computed.
# computing the adjacency matrix of the grid points, with boundaries identified as in the real projective space
from sklearn.utils.graph_shortest_path import graph_shortest_path
# This functions prepares the grid matrix with boundary identification
def make_matrix(rows, cols):
n = rows*cols
M = np.zeros((n,n))
for r in range(rows):
for c in range(cols):
i = r*cols + c
# Two inner diagonals
if c > 0: M[i-1,i] = M[i,i-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# Two outer diagonals
if r > 0: M[i-cols,i] = M[i,i-cols] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# vertical twisted boundary identification
if c == 0: M[n-i-1,i] = M[i,n-i-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# horizontal twisted boundary identification
if r == 0: M[n-i-1,i] = M[i,n-i-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
return M
M = make_matrix(20,20)
# computing the distance matrix of the points over the Klein bottle
rp2 = graph_shortest_path(M)
# Plot of the distance matrix
figure = plt.figure(figsize=(10,10))
plt.imshow(rp2)
plt.title('Reciprocal distance between points over the Klein bottle')
plt.colorbar()
plt.show()
# computing the adjacency matrix of the grid points, with boundaries identified as in the klein bottle
from sklearn.utils.graph_shortest_path import graph_shortest_path
# This functions prepares the grid matrix with boundary identification
def make_matrix(rows, cols):
n = rows*cols
M = np.zeros((n,n))
for r in range(rows):
for c in range(cols):
i = r*cols + c
# Two inner diagonals
if c > 0: M[i-1,i] = M[i,i-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# Two outer diagonals
if r > 0: M[i-cols,i] = M[i,i-cols] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# vertical boundary identification
if c == 0: M[i+cols-1,i] = M[i,i+cols-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
# horizontal twisted boundary identification
if r == 0: M[n-i-1,i] = M[i,n-i-1] = 1 + 0.15*(np.random.rand(1)[0]-0.5)
return M
M = make_matrix(20,20)
# computing the distance matrix of the points over the Klein bottle
klein = graph_shortest_path(M)
# Plot of the distance matrix
figure = plt.figure(figsize=(10,10))
plt.imshow(klein)
plt.title('Reciprocal distance between points over the Klein bottle')
plt.colorbar()
plt.show()
# Saving the results into an array
topological_spaces_mat=np.asarray([rp2, klein])
In the next section we are using giotto to compute the persistent homology groups of the topological spaces we just constructed
# the homology ranks we choose to consider
homologyDimensions = (0, 1 ,2)
persistenceDiagram = hl.VietorisRipsPersistence(metric='precomputed', max_edge_length=10, homology_dimensions=homologyDimensions, n_jobs=-1)
persistenceDiagram.fit(topological_spaces_mat)
# List of all the time-pordered persistent diagrams obtained from the list of correlation matrices
zDiagrams = persistenceDiagram.transform(topological_spaces_mat)
The topological information of the point cloud is synthesised in the persistent diagram. The horizonral axis corresponds to the moment in which an homological generator is born, while the vertical axis corresponds to the moments in which an homological generator dies. The generators of the homology groups (at given rank) are colored differently
# plotting the persistent diagram of the projective space
plot_diagram(zDiagrams[0])
# plotting the persistent diagram of the Klein bottle
plot_diagram(zDiagrams[1])
As you can see from the persistence diagrams, all the betti numbers were found. Some other persistent generators are also appearing, depending on how dense the sampling is and how much noise there is. For example, we see a rahter neat persistent diagram over the Klein bottle (wesee 2 persistent generators for $H^1$ and 1 persistent generator for $H^2$). Notice that all these homology groups are computed over the field $\mathbb{F}_2$.