Skip to content

Algebra of Persistence Modules

Mikael Vejdemo-Johansson edited this page Mar 23, 2024 · 7 revisions

The Barcode Algebra component of the library implements a number of algebraic operations on persistence modules. This page is an attempt to explain what the operations are and how they are implemented. This text will also be reused in a paper MVJ is currently working on.

Finitely Presented Persistence Module

Suppose $G$ and $R$ are two graded sets. We can define the free persistence modules $FG$ and $FR$ as the persistence modules that have interval modules $[|g_i|,\infty)$ or $[|r_i|,\infty)$ for each of the elements in the corresponding basis sets. A linear map $FG \to FR$ determines a persistence module as the cokernel: it is the quotient $FG/img(FR)$. For a matrix representing such a linear map, if it is placed on column-reduced echelon form, that corresponds to a basis change for the relations so that in the new basis, each relation is paired with one specific generator and in the new basis the presentation map is diagonal. The resulting barcode, assuming that $r_i$ is paired with $g_i$, will consist of finite intervals $[|g_i|,|r_i|)$ for each such pair, and infinite intervals $[|g_j|,\infty)$ for any generator not covered by these relations.

The various barcode decomposition theorems (Carlsson-Zomorodian, Crawley-Boevey, and numerous other versions) all basically prove that under sufficiently nice conditions, we can make simultaneous basis changes in all degrees so that these presentation maps are indeed diagonal. While working over a field, we can thus WLOG assume that each presentation can be given as a projection matrix (0-1 matrix with at most one 1 per row and exactly one 1 per column; any relations that don't have any influence on the generators can safely be dropped). However, it may often be useful to track the column reduced echelon form even if it is "merely" triangular and not fully a projection matrix. We shall throughout call such a matrix reduced, and will get ample use of the operation of reducing a matrix. In a persistence setting, we can replace later occurring basis elements with linear combinations with older basis elements -- but we can not replace earlier basis elements using elements that don't yet exist. This requirement gives a direction for any column reduction moves.

In the presence of an ambient module in which both generators and relations embed, there may well be a reason to track all these elements as they are represented in the ambient module. (for an example, think of the chains in a homology computation) These can all be understood as separate maps of free persistence modules and thus represented by separate matrices.

Naming the parts

For the following exposition, we shall assume we have persistence modules $X$, $Y$ presented by graded sets:

  • $G_X = \{x_1,\dots,x_{n_X}\}$,
  • $R_X = \{\xi_1,\dots,\xi_{m_X}\}$,
  • $G_Y = \{y_1,\dots,y_{n_Y}\}$,
  • $R_Y = \{\eta_1,\dots,\eta_{m_Y}\}$

and maps

  • $i_X: R_X\to G_X$
  • $i_Y: R_Y\to G_Y$
  • $f: X\to Y$, which in turn is represented by
  • $\phi: G_X\to G_Y$

For each map we can assume we have access to a matrix that expresses that linear map in the basis sets we have started with.

We sometimes choose in the TDA4j library to represent formal linear sums indexed by a type VectorT and taking coefficients in a type CoefficientT. In this setting, we may think of such formal linear sums as being of type Map[VectorT,CoefficientT] -- with some extra handling to implement arithmetic on these associative maps. Then, the matrices may be (sparsely) represented by lists or sequences of these maps -- or if we have keys of some type enumerating the columns as maps. So a matrix may show up either as a List[Map[VectorT,CoefficientT]] or as a Map[VectorS, Map[VectorT,CoefficientT]]]. We will try to give scala commentary in each section below as well.

Well Definedness of a Map of Persistence Modules

For a map $\phi$ to produce a map of persistence modules, we require:

  • $|x_i| \geq |y_j|$ whenever $\phi_{i,j} \neq 0$. (anything that the bar $x_i$ hits has to already exist when $x_i$ comes into existence)
  • $|\xi_i| \geq |\eta_j|$ whenever $\phi_{i,j} \neq 0$ (assuming the presentation maps are essentially diagonal projection maps as described above -- what we really need is the condition that whenever the images of the $\xi_i$ and $\eta_j$ have a non-zero coefficient in any compositions of maps that sends them into $FG_Y$) -- (for as long as $y_j$ is alive, $x_i$ also has to be alive)

Direct Sum

The simplest of algebraic operations to track and implement is the direct sum of two modules. Observe that for free modules, the direct sum is achieved by taking the disjoint union of the basis sets.

Hence, the direct sum of a pair of persistence modules $X,Y$ is given by the graded sets and maps

  • $G_{X\oplus Y} = G_X\sqcup G_Y = \{x_1,\dots,x_{n_X},y_1,\dots,y_{n_Y}\}$
  • $R_{X\oplus Y} = R_X\sqcup R_Y = \{\xi_1,\dots,\xi_{m_X},\eta_1,\dots,\eta_{m_Y}\}$
  • $i_{X\oplus Y} = \begin{psmallmatrix}i_X & 0 \\\ 0 & i_Y\end{psmallmatrix}$

If we are sure that there is no overlap between the keys addressing any of the basis sets the direct sum can be done by simply taking a union (or concatenating) the corresponding lists and maps. Otherwise, we may want to do something to ensure disjointness of all keys. This would mean something like:

def directSum(iX: Map[BasisX, Map[VectorX, CoefficientT]], iY: Map[BasisY, Map[VectorY, CoefficientT]]): 
  Map[Either[BasisX, BasisY], Map[Either[VectorX,VectorY], CoefficientT]]

Canonical functions $X \to X\oplus Y$, $Y \to X\oplus Y$, $X\oplus Y\to X$, $X\oplus Y\to Y$ can be created by simply tagging (with Right or Left) any map keys to map into the direct sum -- or by ignoring any map keys of one of the sides to map out of the direct sum.

Cokernel

The real work-horse in this part of the library is the cokernel -- especially since by dualizing we can express kernels as cokernels, so by preparing the input we can reuse any code written for cokernels to achieve the kernel computations as well.

The cokernel of a function $f:X\to Y$ is the module $Y/img f$.

Generators of $coker f$ will definitely include the generators of $Y$. Relations of $coker f$ will include the images coming in from $X$ through $\phi$.

However, when $\phi$ maps a basis element $x_i$ to a linear combination of basis elements in $FG_Y$, it may very well happen that some of those basis elements $y_j$ vanish before $x_i$ does. At that point, $x_i$ becomes a linear combination only of the remaining basis elements that the map hits. When that happens, it may be that whatever relation was eliminating something in $Y$ in the cokernel stops eliminating that thing. (TODO: I'm not sure this explanation makes sense entirely) So one such relation could go on to kill off multiple bars in the cokernel.

Translating into operations on free modules, we get that the cokernel has generators $G_Y$, but relations include both $R_Y$ and $G_X$. We start out with a map $i_Y\oplus\phi: R_Y\oplus G_X\to G_Y$, given as a matrix by $\begin{psmallmatrix}i_Y & \phi\end{psmallmatrix}$, but this matrix might no longer be reduced. There might be linear dependencies between relations in $Y$ and generator images in $X$. So in order to establish the actual new persistence module, we would need to reduce the matrix $\begin{psmallmatrix}i_Y & \phi\end{psmallmatrix}$. In the reduced shape, the pivot elements (non-zero entries that are the latest non-zero row in their respective column, and such that no entry in a later occuring column is non-zero after the reduction) give us a pairing of birth times (from the elements in $G_Y$) and death times (from either $R_Y$ or $G_X$ depending on where the pivot is).

Now that we have established $R_{coker f}$ and $G_{coker f}$ from the pivot rows and pivot columns of this reduced matrix, and can take the corresponding submatrix as the presentation map for $coker f$, we would very often also like to have access to a concrete map of modules $Y \to coker f$. In other words, we would like a map $FG_Y\to FG_{coker f}$ that represents the projection onto the quotient.

To get this map, we need to send each basis element in $G_Y$ to its corresponding basis element in $FG_{coker f}$. This would be accomplished by just sending these basis elements by identity maps to the corresponding generators. However, if we (as we likely do) want to have the presentation matrix be a partial diagonal, we will need to make a basis change in $FG_{coker f}$, from the basis given by elements in $G_Y$ to a basis given by the columns in the reduced presentation matrix of the cokernel. For an arbitrary vector in $FG_{coker f}$ we can find its actual image in the cokernel by using the inverse matrix of this presentation matrix and multiplying by that.

Example

For an example, consider the persistence modules $[1,3) \overset{(1 \lambda)}\to [1,2)\oplus[0,3)$.

We have the generators in degrees $|x|=1, |y_1|=1, |y_2|=0$ and relations in degrees $|\xi|=3, |\eta_1|=2, |\eta_2|=3$. Presentation matrices are $(1)$ and $\begin{psmallmatrix}0&1\\1&0\end{psmallmatrix}$. Our function is represented by the matrix $(1 \lambda)$.

The cokernel has generators $y_1, y_2$, and has relations $\eta_1, \eta_2, x$. To get the right result when reducing the presentation matrix, we should ensure that both generators and relations are sorted in degree order: generators are $y_2, y_1$ and relations are $x, \eta_1, \eta_2$. This makes the presentation matrix and its reduction (where as we recall, we are clearing out the youngest entries first):

$$\begin{pmatrix} \lambda & 0 & 1 \\\ 1 & 1 & 0 \end{pmatrix} \Rightarrow \begin{pmatrix} \lambda & -\lambda & 1 \\\ 1 & 0 & 0 \end{pmatrix} \Rightarrow \begin{pmatrix} \lambda & 1 & 0 \\\ 1 & 0 & 0 \end{pmatrix}$$

If we track our operations here, we have a new basis for the relations given by $x$, $\eta_1$, $\eta_2-(1/\lambda)x+(1/\lambda)\eta_1$.

So for the cokernel, we retain the generating set as is, and have a 2-element relations set in the same degrees as the first two elements of the initial relations set: our cokernel relations are in degrees $|\gamma_1| = |x| = 1$ and $|\gamma_2| = |\eta_1| = 2$. We drop the relation in degree $|\eta_2|=3$ since it is reduced to a linear combination of the other two relations.

Now, having cleared everything out in the oldest coefficients at each step, we pair each relationship with the youngest generator in its support that is not already paired with a younger relation. This yields the pairings $\gamma_1\sim y_1$ which produces a barcode $[1,1)$ and $\gamma_2\sim y_2$ which produces a barcode $[0,2)$.

To understand the map from $Y$ to $coker f$, we must figure out how the basis elements $y_1, y_2$ in $G_Y$ are mapped. But we expect to have done a basis change on $FG_{coker f}$, setting the new basis to be the columns of the reduced presentation matrix $\begin{psmallmatrix}\lambda&1 \\ 1&0\end{psmallmatrix}$. We may invert this matrix to get the projection matrix $\begin{psmallmatrix}0&1\\1&-\lambda\end{psmallmatrix}$. But since we drop the first basis element (because its bar vanished), we are left with the bottom row $\begin{psmallmatrix}1&-\lambda\end{psmallmatrix}$. Returning back to the order of generators that the problem started with, this makes the projection map onto the cokernel into the one represented by the matrix $\begin{psmallmatrix}-\lambda&1\end{psmallmatrix}$.

Kernel

There are two approaches to computing the kernel of a map. Either we can first figure out how to compute kernels of maps between free modules, and then set up the computation of the relations for the kernel module as a second such kernel computation. This is the approach in [Skraba - Vejdemo-Johansson] on arXiv.

Alternatively, the kernel dualizes to the cokernel. So we can move to the dual, compute a cokernel there, and translate everything back again.

Moving to the dual means:

  1. Either replacing each barcode endpoint $t$ by $-t$, or replacing the total order of the filtration with its opposite order. Note that in Scala, the opposite order can be quite efficiently accessed by using the way we are using the implicit objects programming techniques to handle orderings to begin with. We make both the kernel and cokernel methods expect an implicit ordering by adding (using ord: Ordering[FiltrationT]) to their call signatures. (this is done implicitly if you declare the FiltrationT type variable as FiltrationT : Ordering) Then we can call cokernel from kernel and give an explicit ordering.
  2. Replacing the matrix representing $\phi$ with its transpose.

Hence, our Scala code for the kernel reduces down to:

  def kernel(
    source: List[PersistenceBar[FiltrationT, AnnotationT]],
    target: List[PersistenceBar[FiltrationT, AnnotationT]],
    matrix: RealMatrix
  )(using
    ord: Ordering[FiltrationT]
  ): List[PersistenceBar[FiltrationT, AnnotationT]] = {
    val dualSource: List[PersistenceBar[FiltrationT, AnnotationT]] =
      target.map(pb => PersistenceBar(pb.dim, pb.upper, pb.lower))
    val dualTarget: List[PersistenceBar[FiltrationT, AnnotationT]] =
      source.map(pb => PersistenceBar(pb.dim, pb.upper, pb.lower))
    val cokernelIntervals =
      cokernel(dualSource, dualTarget, matrix.transpose())(using
        ord = ord.reverse
      )
    cokernelIntervals.map(pb => PersistenceBar(pb.dim, pb.upper, pb.lower))
  }

Notice that most of the work done here is in flipping the positions of lower and upper bounds on the persistence bar data type.

Since the transpose of an inverse is the inverse of the transpose, we can take the transpose of the inverse reduced presentation map to recover the injection map from the kernel into its ambient module.