mtsch / Ripserer.jl

Flexible and efficient persistent homology computation.
https://mtsch.github.io/Ripserer.jl/dev/
MIT License
66 stars 7 forks source link

reconstruct_cycles in more than 2-dimension? #145

Closed chudur-budur closed 2 years ago

chudur-budur commented 3 years ago

Is it possible to implement reconstruct_cycle for more than 2-dimension? Let's say in 3D case, the function will compute all the 2-simplexes that "minimally" surround a 3-dimensional hole.

Also the cycles.jl code has no comment, it's very difficult to understand what is going on, specially in the _find_cycle() function. For example,

best_path = edgetype(g)[]  # what does this line do?
best_sx = first(g.removed) # what does this line do?

Is it possible to add some algorithmic descriptions for the functions in the cycles.jl file? It will be really helpful.

Thanks!! I really appreciate your work.

mtsch commented 3 years ago

Thanks for the kind words, it means a lot!

Sadly the algorithm in question will not work in higher dimensions. It works by first building the 1-skeleton of the simplicial complex at the specified time as a graph. Then, it goes through each edge in the cocycle and connects its endpoints with the shortest path (ignoring the edges of the cocycle). Finally, it selects the one with the minimum weight and returns it (this is what best_sx and best_path are). The problem here is that we don't have a nice analog to a shortest path in higher dimensions.

It might also be worth noting that in general, finding optimal representative cycles in dimensions higher than one is a very hard problem. If I'm not mistaken, the problem of finding the "shortest" representatives is NP-hard. However, I have a hunch that there is an efficient way to get sufficiently nice cycles according to some other metric. This is definitely something I want to explore in the future.

Also, note that this cycle reconstruction feature is very experimental. I agree the code is a bit convoluted and It also has some serious performance issues. I will probably rewrite it later this year when I have more time to work on this stuff again. I'll try to add some comments to the code in the meantime.

By the way, involuted homology (ripserer(...; alg=:involuted)) has no problem working in high dimensions but it computes an arbitrary cycle. Maybe that can be of some use to you?