Flag persistence generators

[ulupo]

This PR add the possibility to retrieve representative simplicies.

• : Compute and retrieve finite representative for dimensions > 0
• : Compute and retrieve essential representative for dimensions > 0
• : Compute and retrieve finite representative for dimension == 0
• : Compute and retrieve essential representative for dimension == 0
• : Add support in the bindings for retrieving representative simplicies
• : Add support in Python to retrieve and output representative simplicies
  • : Document the output regarding representative simplicies
  • : Throw exception if collapser and representative simplicies are enable at the same time
• : Update documentation about representative simplicies
• : Update example basic tutorial to show example of this feature

[ulupo]

If I'm not mistaken, we always return edges that are in the lower diagonal of the distance matrix, i.e. edges [i, j] such that i > j. Is this true consistently? Especially between the case of dimension 0 and of higher dimensions?

[MonkeyBreaker]

Yes, we always return edges [i, j] where we i > j.
Quick question, this order doesn't imply that we are in the upper diagonal instead (i for column and j for row)? Or did I misunderstand something ?

[ulupo]

(i for column and j for row)

We are playing with fire here :) I suggest we stick to i indexing rows and j columns as this is virtually universal in mathematics and Python.

[MonkeyBreaker]

Oh !
Thank you for the correction. Then yes we stick to i for rows and j for columns.

[review-notebook-app]

Check out this pull request on 

See visual diffs & provide feedback on Jupyter Notebooks.

---

Powered by ReviewNB

[ulupo]

@MonkeyBreaker You probably have noticed that there is repeated logic between the 0-dimensional fix and the code for link in the union_find class. An idea would be to make link return what is now called birth_idx instead of being void. Also because I think the check

 if (get_diameter(e) != 0) { 

might be wrong.

[MonkeyBreaker]

@ulupo about returning the index from the link function, it should be doable.
I just need to think to implement it in a way it makes sense when you will read the function.

About the check if (get_diameter(e) != 0) {, I am not sure I understand the issue.
It just ensure that the diameter of the edge is not 0. I am not sure I understand the relation with the vertices, we could maybe have a call ?

[ulupo]

Just a draft suggestion to explain what I mean

[ulupo]

We should add a test that covers sparse matrices with non-zero values along the diagonal. Here we can also test the finite bars in dimension zero and that the birth indices are correct there. Here it is important that the input is valid for a filtration. For example we can create random diagonal entries between 0 and 1, and ensure all off-diagonal entries are greater than 1.

[MonkeyBreaker]

We should add a test that covers sparse matrices with non-zero values along the diagonal. Here we can also test the finite bars in dimension zero and that the birth indices are correct there. Here it is important that the input is valid for a filtration. For example we can create random diagonal entries between 0 and 1, and ensure all off-diagonal entries are greater than 1.

This should be addressed in the latest commit

[ulupo]

Another note: we should avoid using the word "representatives" for what is computed in this PR. We can use "persistence pairs" or "simplex pairs" or "persistence generators"

[MonkeyBreaker]

Well I used "representatives" from the beginning.
But if there is a better name, lets go for it.
I just need to know to update everything.

For me the most important point is the name in the dictionary as output in ripser_parallel. I would like to keep short, 4-5 characters maximum. Any suggestion ?

[ulupo]

I just need to know to update everything.

I took the liberty of renaming ret_representative_simplices to return_generators and renaming the key to 'gens' and rewriting the docs on this part a little.

[ulupo]

Another note: res['gens'][1] and res['gens'][3] should always have length maxdim. Currently, they can be empty lists if there is no content. Instead, we should always return a list of length maxdim with empty ndarrays if necessary in some (or even all) positions.

@MonkeyBreaker I'm guessing this can be fixed at the level of the bindings?

[ulupo]

During work on this PR #29, @MonkeyBreaker and I realized that the current core logic in the matrix reduction algorithm makes it so that we break the following promise to the user: that passing a dense matrix with a few np.inf will lead to exactly the same output as passing a sparse version where the np.inf entries are just missing entries. This was never noticed at the level of barcodes because a "finite" bar with death equal to infinity looks just the same as an infinite bar with the same birth. But flag generators don't lie: we now end up telling the user there are never any generators for essential classes if the matrix is dense, no matter what the matrix is or represents.

Because changing this requires potentially delicate changes to the reduction routine, we decided to discuss what to do about this separately, and open a separate PR if necessary after this one is merged.

[ulupo]

In view of #36 , we need to update the imports in the two notebooks here.

[MonkeyBreaker]

In view of #36 , we need to update the imports in the two notebooks here.

Before doing the merge, we will need to rebase to main and then we will be able to merge, otherwise as you point out, we will have conflicts

[MonkeyBreaker]

Another note: res['gens'][1] and res['gens'][3] should always have length maxdim. Currently, they can be empty lists if there is no content. Instead, we should always return a list of length maxdim with empty ndarrays if necessary in some (or even all) positions.

It seems weird to me, the C++ core has the same logic as for barcodes (barcodes have a maxdim+1 arrays of barcodes because you include dim0). The only difference is that you have maxdim as size of res['gens'][1] and res['gens'][3] , because dimension 0 is not included.

I need to verify if this is not the case

[ulupo]

The only difference is that you have maxdim as size of res['gens'][1] and res['gens'][3] , because dimension 0 is not included.

I don't think this can be the case if maxdim means the same as in Python. Example: if maxdim=1 (compute up to and including dimension 1) then the list has length 1, not 1-1=0.

[MonkeyBreaker]

I don't think this can be the case if maxdim means the same as in Python. Example: if maxdim=1 (compute up to and including dimension 1) then the list has length 1, not 1-1=0.

Well when you compute maxdim=1 you have dimension 0 and dimension 1 => 2-1 = 1, no ?

[ulupo]

Oh I'm sorry @MonkeyBreaker , I thought you said maxdim - 1 (because it used to say that in the docstrings), but you didn't! We said the same thing, i.e. that res['gens'][1] and res['gens'][3] should have length maxdim.

[MonkeyBreaker]

Oh I'm sorry @MonkeyBreaker , I thought you said maxdim - 1 (because it used to say that in the docstrings), but you didn't! We said the same thing, i.e. that res['gens'][1] and res['gens'][3] should have length maxdim.

Does it work as expected then ? Do we have an issue in the documentation ? (In the documentation we are right in that we return maxdim for [1] and [3].

[ulupo]

Does it work as expected then ?

I think I was able to produce examples which returned empty lists when maxdim=1, but can't remember how. Maybe square_dm?

[MonkeyBreaker]

Maybe square_dm?

Just tried with maxdim=1 works as expected