Re: Computing extended persistence in dionysus2

Dmitriy Morozov
 

You can't define extended persistence for an arbitrary filtration. You need a PL function, i.e., you need values on the vertices. Here's a simple example:

import numpy as np
import dionysus as d

w = -1

values = [1,2,3,4]
simplices = [[0], [1], [2], [3], [0,1], [0,3], [1,2], [2,3]]

up_simplices   = [d.Simplex(s,       max(values[v] for v in s)) for s in simplices]
down_simplices = [d.Simplex(s + [w], min(values[v] for v in s)) for s in simplices]
up_simplices.sort()
down_simplices.sort(reverse = True)
f = d.Filtration([d.Simplex([w], -float('inf'))] + up_simplices + down_simplices)

for s in f:
    print(s)

m = d.homology_persistence(f)
dgms = d.init_diagrams(m, f)

for i, dgm in enumerate(dgms):
    print("Dimension", i)
    for pt in dgm:
        print(pt)


This produces the following output. Pair (-inf,inf) is there for technical reasons and should be ignored. Without it, you get the extended persistence diagrams.

<-1> -inf
<0> 1
<1> 2
<2> 3
<3> 4
<0,1> 2
<0,3> 4
<1,2> 3
<2,3> 4
<-1,2,3> 3
<-1,1,2> 2
<-1,0,3> 1
<-1,0,1> 1
<-1,3> 4
<-1,2> 3
<-1,1> 2
<-1,0> 1
Dimension 0
(-inf,inf)
(1,4)
Dimension 1
(4,1)
Dimension 2


On Wed, Oct 3, 2018 at 8:37 PM <cai.507@...> wrote:
Thanks a lot for response. I tried to translate from your example to Dionysus2 and I almost get it but because API for init_diagrams in two versions are different, I fail on the last step.

I try to capture a simple cycle here. I expect the solution should be [3, 0.5] . If I understand coning strategy correctly, The filtration order seems to be right (
<-1> 0
<0> 0.5
<1> 1
<3> 2
<2> 3
<0,1> 1
<0,3> 2
<1,2> 3
<2,3> 3
<-1,2> 0
<-1,3> 0
<-1,1> 0
<-1,0> 0
<-1,2,3> 0
<-1,1,2> 0
<-1,0,1> 0
<-1,0,3> 0)
But the answer it gives is (3,0). I think it's due to I am lacking an analog of eval_ep function in your code. But I looked at the API for init_dionysus in Dionysus2,  and didn't find where to put this eval_ep(value) in.

Hopefully this makes sense to you. Thanks a lot for your time!

Your code (I am having trouble installing Dionysus1 so didn't run the code yet):
# Extended persistence
w = -1
cone = [Simplex([w] + [v for v in s.vertices]) for s in elephant_complex]
cone.append(Simplex([w]))

def projection(points, axis = 1):
def value(v):
return points[v][axis]

return value
value = projection(elephant_points, 1)

def ep_compare(values):
def max_vertex(s):
return max(values(v) for v in s.vertices if v != w)

def min_vertex(s):
return min(values(v) for v in s.vertices if v != w)

def compare(s1, s2):
if s1.dimension() == 0 and w in s1.vertices:
return -1
if s2.dimension() == 0 and w in s2.vertices:
return 1

if s1.dimension() != s2.dimension():
return cmp(s1.dimension(), s2.dimension())

if (w in s1.vertices) ^ (w in s2.vertices): # only one cone simplex
return 1 if w in s1.vertices else -1
elif w in s1.vertices:
return -cmp(min_vertex(s1), min_vertex(s2))
else:
return cmp(max_vertex(s1), max_vertex(s2))

return compare

f = Filtration(elephant_complex + cone)
f.sort(ep_compare(value))

persistence = StaticPersistence(f)
persistence.pair_simplices()

def eval_ep(values):
def eval(s):
if s.dimension() == 0 and w in s.vertices:
return float('inf')
if w in s.vertices:
return min(values(v) for v in s.vertices if v != w)
else:
return max(values(v) for v in s.vertices)
return eval

dgms = init_diagrams(persistence, f, eval_ep(value), lambda i: i)

print show_diagram(dgms)
#print dgms
#print len(dgms)


My code:
import numpy as np
import dionysus as d
w = -1
simplices = [([0], 0.5), ([1], 1), ([2], 3),([3],2), ([0,1], 1), ([0,3], 2), ([1,2],3), ([2,3],3)]

def value(v):
for simplex in simplices:
if len(simplex[0]) == 1 and simplex[0][0]==v:
return simplex[1]

elephant_complex = d.Filtration()
for vertices, time in simplices:
elephant_complex.append(d.Simplex(vertices, time))
cone = [d.Simplex([w] + [v for v in s]) for s in elephant_complex]
cone.append(d.Simplex([w]))

def print_filtration(f):
for simplex in f:
print simplex
def ep_compare(values):
def max_vertex(s):
return max(values(v) for v in s if v != w)

def min_vertex(s):
return min(values(v) for v in s if v != w)

def compare(s1, s2):
if s1.dimension() == 0 and w in s1:
return -1
if s2.dimension() == 0 and w in s2:
return 1

if s1.dimension() != s2.dimension():
return cmp(s1.dimension(), s2.dimension())

if (w in s1) ^ (w in s2): # only one cone simplex
return 1 if w in s1 else -1
elif w in s1:
return -cmp(min_vertex(s1), min_vertex(s2))
else:
return cmp(max_vertex(s1), max_vertex(s2))

return compare

f = d.Filtration(elephant_complex)
for simplex in cone:
f.append(simplex)
f.sort(ep_compare(value))
print_filtration(f)

m = d.homology_persistence(f)
dgms = d.init_diagrams(m, f)

Join dionysus@groups.io to automatically receive all group messages.