Perseus

The Perseus software project for rapid computation of persistent homology.

Download or compile

Download the executable for your operating system. The current version is 4.0 Beta.

Once you have the file for your operating system, see Basic usage below. If you’d prefer to compile the software from source yourself, keep reading.

Download the zipped source to a directory where you have read/write permissions. You can now use any C++ compiler to compile the main file Pers.cpp. The choice of compiler depends mainly on your operating system. Windows users have options such as the open-source MinGW or the complete development environment provided by the somewhat pricey Microsoft Visual Studio. Mac users will probably require a hefty Xcode download and installation.

If using the gcc compiler from the command line, go to the directory with the extracted source files and type:

g++ Pers.cpp -O3 -o perseus

Of course, you can replace perseus above with any executable name of your choice.

Warning for new users: the current version is not compatible with the Clang compiler (the default on new Macs) or with gcc 4.9. For now, new users can either get a genuine gcc for their machines or use the pre-compiled executable. Thanks to Steve and Ken Ferry for pointing out the Clang issue and to Marcio Gameiro for the gcc issue. These will be resolved in the next (Beta 5) version of Perseus.

To report a bug or suggest a new feature, please send me an email. It really helps if the subject line contains the word “Perseus” somewhere!

Basic usage

Perseus computes the persistent homology of many different types of filtered cell complexes after first performing certain homology-preserving Morse-theoretic reductions. In all cases, you prepare the input filtration as a correctly formatted text file (see the formats below) and then read the output persistence intervals, again presented as text files.

Since Perseus is based on discrete Morse theory, it does not rely on the idiosyncrasies of a particular type of complex structure or dimension for its efficiency. That said, some types of complexes are much more common than others in data analysis. For movies and images, it is best to work with cubical data structures. If the data source is a manifold triangulation, the appropriate representation consists of top-cell information on a simplicial complex. Point cloud data is usually handled effectively with Vietoris–Rips complexes built around those points.

To compute the persistent homology of your input, enter the following command-line arguments in this order, either directly in a terminal or in your development environment:

(path to perseus executable) (complex type) (input filename) (output file string)

The (output file string) is optional: it is a prefix for the output files so that they are easy to locate and identify. If you are working on project DeathStar, then making the output string DeathStar produces output files whose names start with DeathStar. To parse the output, see Understanding the output. The sections below explain the acceptable complex types, how to format the corresponding input files, and give simple examples of each. Here is a cheat-sheet of the most common complex types:

cubtop
a dense cubical grid from top-level cube information
scubtop
a sparse cubical grid from top-level cube information
simtop
a simplicial complex from uniform top-level triangulation data
nmfsimtop
a general simplicial complex from non-uniform top-level simplex information
rips
a Vietoris–Rips complex from non-uniform birth point information
brips
a Vietoris–Rips complex from uniform birth point information
distmat
a Vietoris–Rips complex from a distance matrix (corrmat for a correlation matrix)

Perseus is equipped to handle a diverse variety of input complexes, and it is not too difficult to add front-ends for new complex types. There is also no restriction on dimension except system memory: simplices of dimension 70 have a lot of faces and sub-faces! It is important to follow the file formats strictly, so please check the example files in case of confusion.

Warning: while preparing your dataset you will most likely enter “birth times” for various cells. To make sure the software processes your complexes correctly, please avoid using 0 as a birth time. Also do not use floating-point numbers or negative numbers for birth times (except −1, which has a specific purpose). Instead, just use 1, 2, 3, … and so on.

Persistent homology of cubical toplexes

A toplex is a cell complex which can be constructed from a list of top-dimensional cells. So a cubical toplex of dimension d is just a list of d-dimensional cubes; their faces and subfaces are constructed automatically. A simple way to represent a full-dimensional cube is by the coordinate of its lexicographically minimal vertex, which we call its anchor. So a 3D cube anchored at the integer point (a, b, c) is the cube with vertices (a+i, b+j, c+k), where i, j and k take only the values 0 and 1. Perseus can presently compute the persistent homology of the following types of cubical data.

Dense cubical grids

A dense cubical grid of dimension d consists of integer-vertex cubes in d-dimensional Euclidean space. Each cube gets a birth time; if you want standard homology rather than persistence, set the birth times of all cubes that appear to 1. By convention, a birth time of −1 is reserved for cubes that are missing from the grid.

The first line of the file contains the dimension d of the grid. The next d lines give the grid size along each dimension; these must be non-negative. So for a four-dimensional grid of 100 × 30 × 50 × 12 cubes, the first five lines must be 4, 100, 30, 50 and 12, in that order. The remaining lines contain the integer birth times of all the cubes in lexicographical order of anchor points, with preference given to the dimensions in their chosen order.

cubgrid.txt: a 2D 3 × 4 griddownload
2the dimension of the cubical grid
3the number of cubes along the first dimension
4the number of cubes along the second dimension
1birth time of the 2D cube anchored at (0,0)
2birth time of the cube anchored at (1,0)
1birth time of the cube anchored at (2,0)
2birth time of the cube anchored at (0,1)
-1there is no cube anchored at (1,1)
2birth time of the cube anchored at (2,1)
1birth time of the cube anchored at (0,2) … and so on!

In this example the birth time of each 2D cube is the sum of its anchor coordinates, except for the cube anchored at (1,1), which never appears and so gets −1. Note that all indices start at 0. Each lower-dimensional cube automatically inherits the minimum birth time of the higher-dimensional cubes in its coboundary. The keyword is cubtop, short for “cubical toplex”:

./perseus cubtop (path to input file) (output file string)

Since this data structure is stored internally as a contiguous array for constant-time access, each top cube that never shows up (birth time −1) wastes both memory and time. If many cubes never appear, use the sparse format below.

Sparse cubical grids

This format suits a cubical grid where a large fraction of the cubes are missing. If a large hollow sphere or torus is approximated by cubes, we don’t want cubes covering the hollow portions, and the dense format would be very wasteful: an input file almost full of −1s. As a heuristic, say you want the homology of the black region in a 1000 × 1000 black-and-white image. If more than 70 percent of the million pixels are black, the dense format won’t hurt; below 70 percent, you should probably use this sparse format.

The first number is the dimension d of the cubes. After this, each cube is represented by d+1 numbers separated by spaces: the first d are the coordinates of its anchor and the last is its birth time. Avoid −1 as a birth time unless you want that cube to be ignored, in which case there is no need to write it in the file at all.

scubtop.txt: four cubes in dimension 3download
3the dimension of the cubical grid
1 3 4 2the 3D cube anchored at (1,3,4), with birth time 2
2 12 −4 7the 3D cube anchored at (2,12,−4), with birth time 7
0 2 3 1the 3D cube anchored at (0,2,3), with birth time 1

You can of course enter many more cubes. If you enter several cubes with the same anchor and different birth times, only the first is used. In the example, the point (1,3,4) is a face of both the first and the third cube, so it is assigned the lower birth time, 1, from the third. For standard homology rather than persistence, give every cube birth time 1. The keyword is scubtop, for “sparse cubical toplex”:

./perseus scubtop (path to input file) (output file string)

Persistent homology of simplicial toplexes

A simplicial toplex is a list of top-dimensional simplices from which face information is automatically extracted to build a simplicial complex. There are two cases: one where all simplices share a uniform top dimension, and one where they do not. Triangulating a manifold of dimension d always yields uniform top dimension d; a general simplicial complex need not, since it could consist of a 1-simplex sharing a vertex with a 3-simplex.

Representing simplices

First provide the ambient dimension n of the vertices. I strongly recommend indexing all vertices by natural numbers in any order and using 1 as the ambient dimension. This is not necessary, but it is very efficient in time and memory. In general, a point of ambient dimension n is represented by its n floating-point coordinates; with natural-number indexing each point is just a single number.

A simplex of dimension d is represented by its d+1 vertices, so we need n(d+1) numbers per simplex. Here is the 2-simplex with vertices (0.1, 3.4, 7.1), (2.4, −1.6, 0.4) and (13.5, −7.9, 1.0):

0.1  3.4  7.1  2.4  -1.6  0.4  13.5  -7.9  1.0

You can scramble the order of the vertices, but don’t change the order of a single vertex’s coordinates. If instead we index the vertices in lexicographical order, this simplex becomes simply:

1 2 3

Uniform triangulations

The first number in the file is d, the uniform top dimension; the next is n, the number of coordinates per vertex (again, I recommend natural-number labels and n = 1). Each subsequent line contains n(d+1)+1 numbers: the first n(d+1) are the vertices of one simplex, in any order, and the last is an integer birth time. Simplices with birth time −1 are ignored. For ordinary homology, give every simplex birth time 1.

simtop.txt: four 3-simplices with 2 coordinates per vertexdownload
3the dimension of the simplices, so each has 3+1 = 4 vertices
2the number of coordinates per vertex
0 0 0 1 1 0 1 1 1the 3-simplex with vertices (0,0), (0,1), (1,0), (1,1) and birth time 1 … and so on!
./perseus simtop (path to input file) (output file string)

Non-uniform triangulations

Now the toplex has no uniform top dimension. The ambient dimension n of the vertices is still constant, and is the first number in the file (again, I strongly recommend n = 1 with integer labels). Since dimensions vary, each simplex carries its own: before the vertex coordinates, put a number m, the dimension of that simplex, so the software knows that the next n(m+1) numbers are its vertices. The birth time is the last number, as usual. Here, with n = 3, is a 2-simplex with vertices (0.1, 0.3, 7), (−2, 8.1, 0) and (8.9, −3, 9.2) and birth time 5; note the leading 2:

2 0.1 0.3 7 -2 8.1 0 8.9 -3 9.2 5
nmfst.txt: simplices of various dimensionsdownload
1the number of coordinates per vertex
2 1 3 5 1the 2D simplex with vertices 1, 3 and 5; birth time 1
3 1 2 4 6 2the 3D simplex with vertices 1, 2, 4 and 6; birth time 2
6 1 2 3 4 5 6 7 4a 6D simplex with vertices 1 through 7 … and so on

This format allows inconsistencies: you could give a simplex birth time 3 and one of its faces birth time 10. In that case the face’s birth is ruthlessly changed to 3 so that the filtration structure is preserved; no face may be born after its parent simplex. For plain homology, set all birth times to 1. The keyword is nmfsimtop, for “non-manifold simplicial toplex”:

./perseus nmfsimtop (path to input file) (output file string)

Persistent homology of Vietoris–Rips complexes

The Vietoris–Rips complex is completely determined by its 1-skeleton, which can be represented as a symmetric distance matrix whose entries come from pairwise distances in a point cloud, from correlations, or from many other sources. Perseus handles three kinds of data: points with uniform birth times, points with different birth times, and distance matrices.

Points with uniform birth times

This is the most common type of Vietoris–Rips complex. The input is a list of points in Euclidean space, each with an initial radius r. The radius of every point is incremented N times by a universal step size s, giving radii r, r+s, r+2s, …, r+Ns (usually r = 0 for all points, but this is not necessary). Picture balls of growing radius around each vertex: the edge between two vertices is born at the first step between 1 and N at which their balls intersect. If the points are so far apart that their balls never meet, even after growing by Ns, no edge is introduced.

The first line is the ambient dimension d. The next three numbers are, in order, an initial radius multiplier k, the step size s and the number of steps N. The multiplier scales every vertex’s initial radius; almost always you should set k = 1 (with k = 0.5 all initial radii are halved). Then each vertex is a line of d+1 numbers: d coordinates followed by its initial radius, which must be non-negative or the vertex is ignored. For instance, with d = 4, the line 2.3 4.5 -3.7 2.8 0.15 is the vertex (2.3, 4.5, −3.7, 2.8) with initial radius 0.15.

brips.txtdownload
3the ambient dimension, i.e. the number of coordinates per vertex
1 0.01 100radius scaling factor k = 1, step size s = 0.01, number of steps N = 100
1.2 3.4 -0.9 0.5the vertex (1.2, 3.4, −0.9) with radius r = 0.5
2.0 -6.6 4.1 0.3the vertex (2.0, −6.6, 4.1) with radius r = 0.3 … and so on!
./perseus brips (path to input file) (output file string)

Points with different birth times

This type is built from vertex–radius–birth triples. Two vertices get an edge if their balls of the given radii intersect, and the edge gets the maximum of the two birth times. For example, take the triples [(0,0,0); 1; 2] and [(0,1,1); 1; 3]. The distance between the vertices is about 1.414 and the sum of their radii is 2, so there is an edge, born at time 3. Similarly, 2-simplices get the maximum birth time of the edges in their boundaries, and so forth.

The first number n is the embedding dimension of the vertices. Each subsequent line has n floating-point coordinates, a non-negative floating-point radius and an integer birth. Vertices with birth −1 are ignored.

rips.txtdownload
2the ambient dimension, i.e. the number of coordinates per vertex
0.1 0.1 .52 1the vertex (0.1, 0.1) with radius 0.52 and birth 1
-0.2 0.1 .51 3the vertex (−0.2, 0.1) with radius 0.51 and birth 3 … et cetera
./perseus rips (path to input file) (output file string)

Distance and correlation matrices

Sometimes it is infeasible to embed a point sample in Euclidean space. Instead we might have a symmetric distance matrix whose (i, j) entry is the distance between the i-th and j-th points, with no information about their locations. Or we might record the correlations of pairs of events, where the notion of points and distances doesn’t even make sense.

Either way, we can still build a Vietoris–Rips complex. In fact, one of its biggest advantages over its rivals, the Čech and alpha complexes (aside from efficient computability), is that it can be built from a distance matrix rather than point locations. These advantages come at a price: there is no analogue of the nerve lemma for Vietoris–Rips complexes.

First give the size M of the (square) matrix: write 10 for a 10 × 10 matrix. Next give a minimal threshold distance g: any two points less than g apart get an edge born at time zero, so to disconnect all distinct points initially, set g to zero. Then come the step size s and the number of steps N: an edge is born at step n if n is the smallest number between 1 and N for which the threshold ns + g exceeds that pairwise distance; if there is no such n, the edge is never created. Finally, a positive integer C caps the dimension of the complex: no simplices of dimension larger than C are built. This is for computational reasons: if the features of interest are low-dimensional, there is no need to waste time and memory on higher-dimensional complexes. Then enter the M² entries of the matrix; the diagonal entries should be zero.

distmat.txt: M = 3download
3the number of rows/columns of the symmetric distance matrix
0.1 0.2 5 2initial threshold g = 0.1, step size s = 0.2, number of steps N = 5 and dimension cap C = 2
0 0.26 0.4distances from entry 1 to itself, entry 2 and entry 3
0.26 0 2.1distances from entry 2 to entry 1, itself and entry 3
0.4 2.1 0etc.

If you have correlations between event pairs instead of distances, the effective “distance” is 1 − correlation. There is no need to convert your data: Perseus accepts correlation matrices directly and performs the (1 − stuff) calculation itself.

./perseus distmat (path to input distance matrix file) (output file string)
./perseus corrmat (path to input correlation matrix file) (output file string)

Understanding the output: intervals

This section assumes some familiarity with barcodes, persistence diagrams and persistence intervals. Please see the excellent surveys by Carlsson and Ghrist if these terms are unfamiliar. An even more elementary account of simplicial complexes, their geometry and what homology actually measures is in my introductory paper with Sazdanovic, and for a comprehensive (but elementary!) overview of applied topology, consult Rob Ghrist’s book.

Computing persistence diagrams requires working over a field, and all output of Perseus assumes mod-2 coefficients.

Output text files

The main output is a collection of text files of persistence intervals, one per dimension, named using the output file string from the command line. The string can include a whole path: with C:\users\myName\myOutput\result on Windows or /home/myName/myOutput/result on Mac or Unix, the files are created in myOutput with names starting “result”. You need write permission in that directory. With no output string the default is output, and files are created with this prefix in the current working directory. A file with the same name already there will be overwritten!

With the default string, the files are called output_0.txt, output_1.txt, … and so on. How many are created depends on the dimension of the Morse-reduced complex, which is at most the dimension of the input. The reduced version of a 3-dimensional complex might be only 2-dimensional, in which case output_3.txt is not created.

Each file output_n.txt contains the persistence intervals of n-dimensional generators of homology, as two columns of integers: birth time, then death time. A death time of −1 means infinite persistence: the generating cycle persisted past the end of the filtration and was never filled in as a boundary, even once all the cells of the complex had been added.

Finally, output_betti.txt contains the Betti numbers at each step of the filtration. A typical line looks like

12 14 4 7 0

which says that once all cells with birth time at most 12 are included, there are 14 connected components, 4 tunnels, 7 cavities and no higher-dimensional generators of homology. The numbers 14, 4 and 7 are the zeroth, first and second Betti numbers of the 12th subcomplex in the filtration.

Plotting persistence diagrams

To aid visualisation, a simple Matlab script called persdia is bundled with the source code. Call it from the Matlab prompt to plot a Perseus output file as a persistence diagram:

>> persdia('output_1.txt');

You may need to change 'output_1.txt' to the path where your Perseus output files are stored. Here is a sample diagram created by persdia:

A sample persistence diagram: blue birth–death points above a green diagonal, with red diamonds for infinite intervals
The green line is the diagonal and the (birth, death) points are plotted in blue. Points with infinite persistence (death −1 in the output file) are plotted as red diamonds above their birth times.

License and citing

Perseus 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. 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.

If Perseus was useful in your research, please cite it! The Morse-theoretic simplification engine of Perseus comes from this paper:

Konstantin Mischaikow and Vidit Nanda. Morse Theory for Filtrations and Efficient Computation of Persistent Homology. Discrete & Computational Geometry, 50(2): 330–353, September 2013.

To cite the software itself, please use:

Vidit Nanda. Perseus, the Persistent Homology Software. https://people.maths.ox.ac.uk/nanda/perseus/, accessed DD/MM/YY.

And finally, you might try to locate yourself on this map of the parts of the world where Perseus was recently downloaded:

World map highlighting places where Perseus has been downloaded