domain: info math

text:
Peter Bubenik
email: peter.bubenik@gmail.com

L D Pawe

Lotko

A PERSISTENCE LANDSCAPES TOOLBOX FOR TOPOLOGICAL STATISTICS

Keywords: topological data analysis, persistent homology, statistical topology, topological machine learning, intrinsic dimension

Topological data analysis provides a multiscale description of the geometry and topology of quantitative data. The persistence landscape is a topological summary that can be easily combined with tools from statistics and machine learning. We give efficient algorithms for calculating persistence landscapes, their averages, and distances between such averages. We discuss an implementation of these algorithms and some related procedures. These are intended to facilitate the combination of statistics and machine learning with topological data analysis. We present an experiment showing that the low-dimensional persistence landscapes of points sampled from spheres (and boxes) of varying dimensions differ.

Introduction

We provide some algorithms and computational tools for statistical topological data analysis. In particular, we give algorithms for calculating the persistence landscape, a functional summary of persistence modules. We also give algorithms for calculating the averages of such summaries, and for calculating distances between such averages. These tools also provide an alternative computational approach for calculating distances between topological summaries that may be useful when other methods are computational prohibitive. In addition, we specify an implementation of these algorithms and some related tools that we have made publicly available.

We are motivated by topological data analysis [START_REF] Ghrist | Barcodes: the persistent topology of data[END_REF][START_REF] Carlsson | Topology and data[END_REF]. Its main tool, persistent homology provides a multiscale description of the topology of the data of interest, called either a barcode or a persistence diagram. Unfortunately this summary is difficult to work with from the point of view of statistics and machine learning. For example, it is not feasible to calculate averages. For these purposes, it is convenient to replace these summaries with a linear summary, that is, a finite-or infinite-dimensional vector. In a linear space it is easy to calculate averages. One such vector which does not lose any information is the functional summary called the persistence landscape [START_REF] Bubenik | Statistical topological data analysis using persistence landscapes[END_REF]. Since this summary may be thought of as lying in a Hilbert space, in the language of machine learning, it is a feature map. There is an associated kernel [START_REF] Reininghaus | A stable multi-scale kernel for topological machine learning[END_REF] to which standard machine learning tools may be applied.

1.1. Background. In the simplest computational setting for topological data analysis, the data of interest is encoded in a finite filtered complex, (1)

K 0 ⊂ K 1 ⊂ . . . ⊂ K n .
This is a filtration of the complex K = K n and it is sometime convenient to add K -1 = ∅.

Persistent homology [START_REF] Edelsbrunner | Topological persistence and simplification[END_REF][START_REF] Zomorodian | Computing persistent homology[END_REF] gives a multiscale representation of the topology of this complex. To be precise, one applies homology in some degree with coefficients in some field to [START_REF] Adcock | The ring of algebraic functions on persistence bar codes[END_REF] to obtain a sequence of finite-dimensional vector spaces and linear maps,

(2) H(K 0 ) → H(K 1 ) → . . . → H(K n ), called a persistence module. It turns out that the persistence module can be completely described by a finite sequence of pairs {(b i , d i )}, with b i < d i . For each such pair (b i , d i ) there is a choice of a nonzero homology class α i ∈ H(K bi ) that is not in the image of H(K bi-1 ) and whose image is nonzero in H(K di-1 ) but is zero in H(K di ). One sometimes says that α i is born at b i and dies at d i . Furthermore, the homology classes {α i } and their nonzero images under the maps in [START_REF] Pankaj | Vertical decomposition of shallow levels in 3-dimensional arrangements and its applications[END_REF] give a basis for the vector spaces in [START_REF] Pankaj | Vertical decomposition of shallow levels in 3-dimensional arrangements and its applications[END_REF]. Considering these pairs as points in the plane, one obtains the persistence diagram. Considering them as intervals [b i , d i ) on obtains the barcode. We will often refer to them as birth-death pairs. In the simple setting of (1), we have b i , d i ∈ {0, 1, . . . , n}.

However we can generalize to b i , d i ∈ R by associated a corresponding increasing sequence of real numbers with [START_REF] Adcock | The ring of algebraic functions on persistence bar codes[END_REF]. This summary is stable [START_REF] Cohen-Steiner | Stability of persistence diagrams[END_REF][START_REF] Cohen-Steiner | Lipschitz functions have lp-stable persistence[END_REF][START_REF] Chazal | The structure and stability of persistence modules[END_REF] in that small perturbations of the data will lead to small perturbations of these pairs, under suitable choices of distance. Successful applications of topological data analysis include breast cancer data [START_REF] Nicolau | Topology based data analysis identifies a subgroup of breast cancers with a unique mutational profile and excellent survival[END_REF], sensor networks [START_REF] Silva | Coverage in sensor networks via persistent homology[END_REF], orthodontic data [START_REF] Gamble | Exploring uses of persistent homology for statistical analysis of landmark-based shape data[END_REF], signal analysis [START_REF] Perea | Sliding windows and persistence: An application of topological methods to signal analysis[END_REF], target tracking [4], and brain artery data [START_REF] Bendich | Persistent homology analysis of brain artery trees[END_REF]. Now let us define the persistence landscape [START_REF] Bubenik | Statistical topological data analysis using persistence landscapes[END_REF]. First, for a birth-death pair (b, d), let us define the piecewise linear function f (b,d) : R → [0, ∞].

(3)

f (b,d) =      0 if x ∈ (b, d) x -b if x ∈ (b, b+d 2 ] -x + d if x ∈ ( b+d 2 , d)
The persistence landscape of the birth-death pairs {(b i , d i )} n i=1 is the sequence of functions λ k : R → [0, ∞], k = 1, 2, 3, . . . where λ k (x) is the k-th largest value of {f (bi,di) (x)} n i=1 . We set λ k (x) = 0 if the k-th largest value does not exist; so λ k = 0 for k > n. Equivalently, the persistence landscape is a function λ : N × R → [0, ∞], where λ(k, t) = λ k (t). In this definition we have assumed that b and d are finite. In the appendix we show that this definition extends to the cases where b and/or d are infinite.

Given a set of persistence landscapes, λ (1) , . . . , λ (N ) , their average, λ, is defined pointwise, λk (t) = 1 N N i=1 λ (i) k (t). Distances between persistence landscapes and between average persistence landscapes can be given using the L ∞ norm,

λ -λ ∞ = sup k,t |λ k (t) -λ k (t)| , or the L p norm, for 1 ≤ p < ∞, λ -λ p = ∞ k=1 |λ k (t) -λ k (t)| p dt 1 p
.

In [START_REF] Bubenik | Statistical topological data analysis using persistence landscapes[END_REF] it is shown that the persistence landscape is stable with respect to the L p distance for 1 ≤ p ≤ ∞. That is, under suitable hypotheses, sufficiently small perturbations of a function under the supremum norm lead to small changes of the persistence landscape of the persistent homology of the sublevel sets of that function under the L p norm.

In addition to the persistence landscape, other functional summaries of the persistence module can also be easily averaged and used for statistics and machine learning. The simplest of these is the ranks of the maps in the persistence module (called the persistent Betti number function or the rank function). To help keep distance finite, one can smooth the persistence diagram [START_REF] Donatini | Size functions for signature recognition[END_REF][START_REF] Ferri | Point selection: a new comparison scheme for size functions (with an application to monogram recognition)[END_REF][START_REF] Chepushtanova | Persistence images: An alternative persistent homology representation[END_REF] or integrate with respect to a weight function [START_REF] Robins | Principal component analysis of persistent homology rank functions with case studies of spatial point patterns, sphere packing and colloids[END_REF]. The most sophisticated such smoothing is given in [START_REF] Reininghaus | A stable multi-scale kernel for topological machine learning[END_REF]. Variants of the persistence landscape such as silhouettes [START_REF] Chazal | On the bootstrap for persistence diagrams and landscapes[END_REF] also work. The representation of persistence diagrams as complex polynomials [START_REF] Fabio | Comparing persistence diagrams through complex vectors[END_REF] has been used for shape classification. In [START_REF] Carrire | Stable topological signatures for points on 3d shapes[END_REF] the authors construct stable (with respect to the Bottleneck distance) feature vectors based on the persistence diagram. Those vectors are used later to compare points in 3d shapes. In addition, one can use the lengths of the N longest bars [START_REF] Bendich | Persistent homology analysis of brain artery trees[END_REF] or algebraic functions on barcodes [START_REF] Adcock | The ring of algebraic functions on persistence bar codes[END_REF].

1.2. Relation to other software. The main aim of most persistent homology software is to start with a finite filtered complex (or perhaps something used to generate such as a complex, such as a finite set of points in Euclidean space) and to produce a set of birth-death pairs. Examples include JavaPlex [START_REF] Sexton | PLEX libdary[END_REF], Dionysus [START_REF] Morozov | The Dionysus software project[END_REF], Perseus [START_REF] Nanda | The Perseus software project[END_REF], PHAT [START_REF] Bauer | PHAT (Persistent Homology Algorithm Toolbox[END_REF], and GUDHI [START_REF] Maria | Simplicial Complexes and Persistent Homology Packages[END_REF]. These birth-death pairs are the input for our software.

Other recent software developed concurrently to ours is the R package TDA [START_REF] Fasy | Introduction to the r package tda[END_REF]. This package provides R users with tools for topological data analysis. For example, it provides an interface to GUDHI, Dionysus and PHAT. It calculates persistence landscapes using grids and also calculates confidence intervals. Future versions of TDA may include an interface to our code.

The software described in this paper performs exact computations of persistence landscapes. It also implements grid-based computations, such as the ones provided in the TDA [START_REF] Fasy | Introduction to the r package tda[END_REF] package. Since the grid-based computations are straightforward, they are not discussed in this paper. One may switch from the (default) exact computations to the grid-based ones, by making a small change in the self-explaining config file in the library's main folder. Pros and cons of grid and exact computations of landscapes along with some description of the limitations of both approaches are provided in Section 6.

1.3. Prior work and related work. It has been shown [START_REF] Cohen-Steiner | Lipschitz functions have lp-stable persistence[END_REF] that the set of persistence diagrams with the Wasserstein distance is a complete and separable metric space, and thus provides a suitable setting for probability and statistics. Unfortunately the Fréchet mean is not necessarily unique. For a slightly adjusted metric, there is an algorithm [START_REF] Turner | Fréchet means for distributions of persistence diagrams[END_REF] that converges to an element of the Fréchet mean set, though it does not have good computational properties. The discontinuity of this procedure can be remedied by using a probabilistic approach [START_REF] Munch | Probabilistic fréchet means and statistics on vineyards[END_REF].

Persistence diagrams can be used for statistical inference. Hypothesis testing using persistence diagrams has been considered for brain MRI data in [START_REF] Chung | Persistence diagrams in cortical surface data[END_REF] and more abstractly in [START_REF] Robinson | Hypothesis testing for topological data analysis[END_REF]. Furthermore, confidence sets for persistence diagrams have been obtained in [START_REF] Fasy | Confidence sets for persistence diagrams[END_REF].

The persistence landscape allows the use of more statistical machinery. The bootstrap has been applied to obtain confidence bands for the persistence landscape [START_REF] Chazal | On the bootstrap for persistence diagrams and landscapes[END_REF] and the average persistence landscape of subsamples has been studied [START_REF] Chazal | Subsampling methods for persistent homology[END_REF]. The persistence landscape has been used to study protein binding [START_REF] Nikolic | Using cycles in high dimensional data to analyze protein binding[END_REF] and as a kernel for topological machine learning and compared to the recent multi-scale kernel for persistence diagrams [START_REF] Reininghaus | A stable multi-scale kernel for topological machine learning[END_REF].

1.4. Our work. We present an algorithm for calculating the persistence landscape corresponding to n birth-death pairs in time O(n 2 ) and show that this time complexity is optimal. Given N persistence landscapes, each obtained from n birth-death pairs, we calculate their average persistence landscape in time O(n 2 N log N ). We show how to calculate the L p distances for 1 ≤ p ≤ ∞ between two such averages in time O(n 2 N log N ).

If we are willing to slightly perturb our birth and death times, then we give algorithms that are much faster for large n or large N . To be precise, instead of associating an arbitrary increasing sequence of real numbers with each of our filtered simplicial complex, we round these numbers so that they lie on an equally-spaced grid of size m. By the stability theorem of [START_REF] Cohen-Steiner | Stability of persistence diagrams[END_REF], this will perturb the resulting persistence diagram by at most δ 2 in the bottleneck distance, where δ is the spacing of the grid. Under this assumption, we give an algorithm for calculating the persistence landscape in time O(mn log n) and calculating the average landscape in time O(mnN ). We can also calculate the distance between two such averages in time O(mnN ).

We describe an implementation of these algorithms that we have made publicly available. In addition to the above algorithms, we have also implemented a number of procedures that we hope will be helpful to practitioners interested in using these methods for topological data analysis. For example, we have procedures for plotting persistence landscapes and their averages. Given a number of classes of birth-death pairs, we allow the user to calculate the distance matrix of the corresponding average persistence landscapes, and also to perform pairwise permutation tests for these classes, using the L p distance, with 1 ≤ p ≤ ∞, between their respective average persistence landscapes. In addition, we provide two nearest-neighbor classifiers, one using the persistence landscape in a single degree, and the other using multiple degrees. Also, one can compute the inner product of landscapes with the presented software, thus allowing them to be used with kernel methods.

In Section 2 we describe the input and output for our algorithms. In Section 3 we give our main algorithms and calculate their time complexities. We also show why calculating the persistence landscape is not the same as the n-th envelope problem. In Section 4 we give a few simple numerical experiments demonstrating our implementation of our algorithms. In Section 5 we describe our implementation of our algorithms and some related procedures.

Data structures

In this section, we describe the inputs and outputs of our main algorithms.

2.1. Input. The initial input to our algorithms consists of a list of n pairs of numbers (b, d) with b < d. Each of these pairs represents the birth and death times of a persistent homology class. Thinking of these pairs as points we obtain a persistence diagram [START_REF] Edelsbrunner | Computational Topology[END_REF], and considering them to be intervals we obtain a barcode [START_REF] Ghrist | Barcodes: the persistent topology of data[END_REF]. We will give two algorithms for calculating persistence landscapes from these birth-death pairs.

In Algorithm 1, for simplicity, we assume that b and d are finite. In the appendix, we extend this algorithm to Algorithm 4 which also accommodates infinite intervals. Reduction to the finite case can be achieved by removing or truncating infinite intervals or by using extended persistence [START_REF] Cohen-Steiner | Extending persistence using poincare and lefschetz duality[END_REF], where extended persistence is used to obtain a persistence module which is eventually zero [START_REF] Bubenik | Categorification of persistent homology[END_REF].

In our implementation, the user is asked to define a number i (which can be set to the maximal representable double number) at which infinite intervals will be truncated. This choice will depend on the application and/or the persistent homology calculation. For example, the persistent homology calculation of a filtered Vietoris-Rips complex is often truncated at some filtration value, because of the exponential growth in size of the complex. It is then sensible to truncate infinite intervals at this maximum filtration value.

In Algorithm 2, we assume that each b and d is an element of a finite, evenly-spaced grid, a, a + d, a + 2d, . . . , a + md. Such birth-death pairs are often the output of Perseus [START_REF] Mischaikow | Morse theory for filtrations and efficient computation of persistent homology[END_REF][START_REF] Nanda | The Perseus software project[END_REF] or Plex [START_REF] Sexton | PLEX libdary[END_REF]. In fact, by rescaling, we assume without loss of generality that values of b and d are elements of {0, 2, 4, . . . , 2m}.

2.2.

Output. In this section, we describe our encoding of persistence landscapes and linear combinations of persistence landscapes.

As defined in Section 1, a persistence landscape is a function λ :

N × R → [0, ∞],
or equivalently, a sequence of functions λ k : R → [0, ∞] where k ≥ 1. For every fixed k, λ k is a piecewise-linear function.

Since the input consists of n birth-death pairs, λ k = 0 for k greater than some fixed K ≤ n. We will represent λ k by a vector L k of the points (x, λ k (x)) such that λ k is not differentiable in x, which is sorted to have increasing values of x. For clarity, we include in L k the points (-∞, 0) and (∞, 0). The projection of L k onto its first coordinate is the vector of critical numbers and the projection on the second coordinate is the vector of critical values. We refer to the elements of L k as critical points. Clearly, λ k can be recovered from L k by linearly interpolating consecutive points. This is therefore an exact representation of a persistence landscape. Therefore from now on these two objects will be used interchangeably. Note that when a landscape is represented on a discrete grid, this property do not hold in general.

Let P denote the total number of critical points in the L k , not including the points (±∞, 0), which are not strictly necessary. In Section 3.2 we show that P = O(n 2 ).

In Algorithm 2, the input numbers are in the set {0, 2, 4, . . . , 2m}. It follows that the critical numbers and critical values are elements of {0, 1, 2, . . . , 2m}. We represent the persistent landscape by a two-dimensional array V satisfying V [k][i] = λ k (i), with i ∈ {0, 1, 2, . . . , 2m} and k ∈ {1, . . . , K}, where K is the largest k such that λ k is not identically equal to 0. Note that the persistence landscape can be obtained from V by linear interpolation. It is also fruitful to consider V as a vector of size K(2m + 1).

Notice that a linear combination of persistence landscapes is also a sequence of piecewise linear functions. So we encode it in the same way as we do a persistence landscape.

Remark 1. We remark that K is the largest rank of the linear maps in (2).

Algorithms

In this section, we describe our main algorithms, which compute persistence landscapes, linear combinations of such persistence landscapes, and distances between such linear combinations.

3.1. Persistence landscape. In this section, we present two algorithms to construct the persistence landscape from a list of birth-death pairs. For simplicity, Algorithm 1 assumes the input consists of finite numbers. For a variation without this assumption and with the same complexity, see Algorithm 4 in the appendix. The computational complexity is O(n log(n) + nK), where n denotes the number of input pairs and K the number of nonzero landscapes. Since K ≤ n, this algorithm is O(n 2 ). If we do not need all of the persistence landscape, then there are faster variations, described below. Algorithm 2 assumes that the input coordinates are elements of a finite, evenly-spaced grid of size m + 1. Its computational complexity is O(mn

+ mK log K) which is O(mn log n).
An example of the steps in Algorithm 1 is given in Figure 1.

In Algorithm 1, we first sort the input list, A, in a way that we will be able to pass through A exactly once in order to construct each list L k . This takes time O(n log n). We denote the last pair in L k by L k . last. Each iteration of the outer while loop constructs one of the lists L k . Denote the number of these by K. The length of A is initially n. At the start of the outer while loop, the length of A is decreased by one. The inner while loop does not increase the length of A. So K ≤ n. In each non-terminal iteration of the inner while loop, the position of p in the list A advances, so it repeats at most n times. Thus the algorithm terminates. Furthermore, in each iteration of the outer while loop, each pair in A is only considered once, so it takes time O(n). Therefore, Algorithm 1 takes time O(n log n + Kn).

The following faster variations of Algorithm 1 might be of interest for some applications. If we only want L 1 , . . . , L log n , then the algorithm takes time O(n log n). If we only want L 1 , . . . , L κ , for some fixed κ, then the algorithm also takes O(n log n). If in addition A is already sorted, then the algorithm takes O(n).

In Algorithm 2, we add each birth-death pair's contribution to the persistence landscape to an array of lists W . This takes time O(mn). Next we sort each of the lists, which takes time O(mK log K). Finally we copy this data to V , which takes time O(mK). So the algorithm has time complexity O(mn + mK log K). Since K ≤ n, this is O(mn log n).

3.2.

Persistence landscapes and the n-th envelope problem. The problem discussed in this paper, may at first glance may look like the problem of finding envelopes of a set of line segments, which is well known in computational geometry. In this section, we discuss the finding-envelopes problem from computational geometry and show why the problem of computing the persistence landscape is different. Also we will give a worst case estimate of the spatial complexity of the Algorithm 1: Compute the persistence landscape.

Input: A = {(bi, di)} n i=1 -a list of birth-death pairs, -∞ < bi < di < ∞. Output: {L k } -the
if b > d then Add (d, 0) to L k ; end if b ≥ d then Add (b , 0) to L k ; else Add ( b +d 2 , d-b 2 ) to L k ; Push (b , d) into A in order, starting at p; Let p point to the next term; end Add ( b +d 2 , d -b 2 ) to L k ; (b, d) ← (b , d ); end end ++k; end Return {L k };
persistence landscape which will later guarantee the optimality of the algorithms used in this paper.

The upper envelope of a set of line segments {L i } k i=1 is defined as those parts of the line segments which are visible from the point (0, +∞). Equivalently if we consider L i to be a piecewise-linear functions set to -∞ outside the line segments, the upper envelope in the point x ∈ R is defined as max i∈{1,...,k} L i (x). Due to the practical importance of this problem, there are many efficient algorithms available to compute upper envelopes [START_REF] Hershberger | Finding the upper envelope of n line segments in o(n log n) time[END_REF] practically in linear time. Let us have a family of persistence intervals {(a i , b i )} n i=1 . For every interval let us define the line segment L 2i-1 starting at (a i , 0) and ending at ( ai+bi 2 , bi-ai 2 ) and the line segment L 2i starting at ( ai+bi 2 , bi-ai 2 ) and ending at (b i , 0). The first landscape λ 1 is the upper envelope of the family {L i } 2n i=1 . However, from the classical upper envelope computations we are not able to get the λ k 's for k > 1. There are algorithms available in which successive upper envelopes are computed, as presented in [START_REF] Hershberger | Upper envelope onion peeling[END_REF]. However in this case, once part of a line segment L i belongs to the k-th envelope, then the whole of L i is removed from the search of (k + 1)-st envelope, which is not what is required in the calculation of the persistence landscape. Alternatively, one can find all the intersections between the mentioned lines. However, in this case, we will get O(n 2 ) lines segment in the worst case.

Even though standard results on the n-th envelope problem do not apply directly, Michael Kerber has recently pointed out to us that an adaptation of a line sweep algorithm can indeed be applied to calculate persistence landscapes. The time complexity of this adaptation may be a bit worse than the one presented in this paper, but there are examples for which it is faster. This approach will be explored and compared to one presented here in a subsequent paper.

To end this section, let us determine the complexity of the persistence landscape structure. Clearly the size of the structure and therefore the computational time required to construct the structure is bounded from below by P , the number of critical points of persistence landscape. This number can be quadratic in n in the worst case as shown in the Figure 2. This analysis shows that the algorithms presented in this section are optimal in the worst case. In fact, it is easy to see that if n is the number of intervals given by the birth-death pairs, and p is the number of pairwise non-empty intersections of these intervals, then P = 3n + 2p. Since p ≤ n 2 , we have P ≤ n 2 + 2n.

3.3.

Averages and linear combinations. In this section, we present algorithms for calculating linear combinations of persistence landscapes.

Algorithm 2: Compute the persistence landscape using a grid.

Input: {(bi, di)} n i=1 -a list of pairs bi < di which are elements of 0, 2, 4, . . . , 2m; Output: V -the persistence landscape, a two-dimensional array of size K × (2m + 1); Initialize W , an array of size 2m + 1 of empty lists of integers; Now suppose that we have persistence landscapes λ 1 , . . . , λ N and we wish to calculate the linear combination f = N j=1 a j λ j , where a j ∈ R. Important special cases are the average of a list of persistence landscapes, λ = N j=1 1 N λ j , and the difference between the averages of two groups of persistence landscapes, λ -λ = N j=1

for i = 1 to n do for j = 1 to d i -b i 2 do Append j to W [bi + j]; end for j = 1 to d i -b i 2 -1 do Append j to W [di -j]; end end for i = 0 to 2m do Sort W [i] in decreasing order; end K ← max i∈{0,...,2m} length of W [i]; Initialize V as a K × (2m + 1) zero matrix; for i = 0 to 2m do for k = 0 to length of W [i] do V[k][i] = W[k][i]; end end
1 N λ j + N j=1 -1 N λ j . Let f k (t) = f (k, t). Then f k = N j=1 a j λ j k .
First we give a naive algorithm for calculating a representation of

f k from the representations (X 1 k , Y 1 k ), . . . , (X N k , Y N k ) of λ 1 k , . . . , λ N k . First we sort
Algorithm 3: Linear combination of persistence landscapes.

Input: (X 1 k , Y 1 k ), . . . , (X N k , Y N k )
, for some fixed k, and a1, . . . , aN ; Output: (X k , Y k ) ; Merge the sorted lists X j k , removing duplicates. Call this vector X k ;

for j = 1 to N do Calculate Ȳj k = λ j k (X k ) by linear interpolation; end Y k ← N j=1 aj Ȳj k ; return (X k , Y k )
the union of the elements of X 1 k , . . . , X N k , removing repetitions. Call this vector X k . For each 1 ≤ j ≤ N , define Ȳj k = λ j k (X k ). Now we represent λ j k by X k and Ȳj k . Again, λ j k can be obtained from X k and Ȳj k by linear interpolation. By definition,

f k (X k ) = N j=1 a j λ j k (X k ) = N j=1 a j Ȳj k .
Also, f k can be recovered from X k and f k (X k ) by linear interpolation. See Algorithm 3. In summary, vector space operations on λ 1 k , . . . , λ N k are obtained from vector space operations on Ȳ1 k , . . . , ȲN k . Let us consider the time complexity of this algorithm. Let n be the maximum number of birthdeath pairs used to construct each of λ 1 , . . . , λ N . Let P k be the sum of the number of the critical points of the λ 1 k , . . . , λ N k . Note that P k = O(N n). Then f k has at most P k critical points, and constructing X k takes O(P k ). Since the length of X k may be at most P k , calculating each Ȳj

k takes O(P k ) and calculating Y k = f x (X k ) takes O(N P k ). So Algorithm 3 has time complexity O(nN 2 ).
Now if we repeat Algorithm 3 for all k we obtain a linear combination of the full persistence landscape. Let P be the sum of the number of critical points of λ 1 , . . . , λ N . Then constructing

{L k } has time complexity O(N P ) = O(n 2 N 2 ).
The complexity of the naive algorithm can be improved to O(n 2 N log N ) by merging landscapes in a binary tree fashion (this is sometimes referred to as divide-and-conquer). Here we describe this for the case of calculating the average landscape 1 . Suppose a collection of persistence landscapes

(X 1 k , Y 1 k ), . . . , (X N k , Y N k ) is given. Let us assume that N is even. Then this collection is transformed to a new collection (M 1 X 1 k , M 1 Y 1 k ), . . . , (M 1 X N 2 k , M 1 Y N 2 k ) where (M 1 X i k , M 1 Y i k ) = (X i k , Y i k ) + (X i+1 k , Y i+1 k ) for i ∈ {1, . . . , N 2 }. In the case that N is odd, additionally set (M 1 X N +1 2 k , M 1 Y N +1 2 k ) = (X N k , Y N k )
. This merging procedure is repeated for the obtained sequence until the sequence contains only one element (X f k , Y f k ). The average persistence landscape is then 1 N (X f k , Y f k ). Note that the number of landscapes at level i is of the order N 2 i . The complexity of each of these is 2 i n 2 . Therefore, the cost of merging neighboring persistence landscapes at level i is 2 i n 2 N 2 i = n 2 N , which does not depend on i. Number of levels is log 2 (N ) , and therefore the total cost of computing average persistence landscape is O(n 2 N log(N )).

For the case where all of the birth-death pairs have endpoints on an evenly-spaced grid of size m + 1, we can obtain a linear combination of the persistence landscapes by simply taking a linear combination of the N corresponding vectors of size K × (2m + 1), which has time complexity O(KmN ) which is O(mnN ).

3.4. Distances. In this section, we compute the L p and L ∞ distances between two linear combinations of persistence landscapes L = {L k } and L = {L k }. Let P be the maximum number of critical points of L and L . Let {(X k , Y k )} be the representation of the difference between these two persistence landscapes as described in Section 3.3. 1 In the general case of computations of weighted sums of landscapes, when merging in the first level we compute weighted sums (M

1 X i k , M 1 Y i k ) = a i (X i k , Y i k ) + a i+1 (X i+1 k , Y i+1 k ) instead of unweighted ones (M 1 X i k , M 1 Y i k ) = (X i k , Y i k ) + (X i+1 k , Y i+1 k )
. Also, at the end of computations we do not multiply the last landscape by 1 N .

The L ∞ distance between L and L is

L -L ∞ = max k |Y k | .
This calculation has time complexity O(P ). The L p distance between L and L is given by the formula:

L -L p = K k=1 L k -L k p p 1 p
The norm L k -L k p p can be computed from (X k , Y k ) by summing integrals over intervals given by consecutive elements of X k . These have the form Now assume that we start with two persistence diagrams, each of which has at most n birthdeath pairs. Since P = O(n 2 ) we can calculate the L ∞ and L p distances between the corresponding persistence landscapes in O(n 2 ).

For the case where the points in the birth-death pairs lie on an evenly-spaced grid of size m, the calculation of distance between linear combinations of persistence landscapes is O(Km), where K be the maximum k for which either {L k } or {L k } is nontrivial. So starting with two persistence diagrams, each of which has at most n birth-death pairs which lie on a grid of size m, we can calculate the L ∞ and L p distances between their persistence landscapes in O(Km) which is O(mn).

We can combine the results of this section with those of the previous section to calculate the L ∞ and L p distances between f and g, two linear combinations of at most N persistence landscapes, each of which in obtained from at most n birth-death pairs. An important special case of this is the calculation of the distance between two average persistence landscapes. Let P be the total number of critical points in the two sets of persistence landscapes. Then the representation of f -g, {(X k , Y k )}, has at most P critical points, and can be constructed in time O(n 2 N log N ).

From this, the distance between f and g can be calculated in O(P ) = O(n 2 N ) so the total computation time is O(n 2 N log N ). For the case were the endpoints of the intervals lie on a grid of size m + 1, constructing f and g is O(mnN ), and calculating the distance between f and g is O(mn), so the total computation has time O(mnN ).

Experiments

In this section, we present the results of some experiments used to test our implementation of our algorithms.

4.1.

Points sampled from S d . Our first experiment is motivated by the following questions. Suppose we are given a set of points in R D that lie on a lower dimensional sphere S d . How well can we determine d using persistent homology?

For d ∈ {2, 3, . . . , 10}, we sampled 100 points from S d using the uniform distribution. This was done by sampling 100 points from the (d + 1)-dimensional Gaussian distribution and then projecting those points to S d . For an example of such normalization in topological data analysis, see [START_REF] Carlsson | On the local behavior of spaces of natural images[END_REF]. This was repeated 1000 times. To refute the charge that our detected differences are only due to the increase of the average distance between points in higher dimensions, we rescale the sphere on which the points has been projected so that the average distance between points is one.

For each point cloud generated as described above we computed the persistent homology of the corresponding Vietoris-Rips complex. The parameter values for the radius ranged from 0 to a radius for which all inessential 0, 1, and 2 dimensional cycles are killed, which was 0.7 for this range of dimensions.

Remark 2. In this example and all subsequent examples, the filtration values have been rescaled to the range from 0 to 100. For simplicity, we have chosen to leave all subsequent calculation in this new scale. If desired, it is easy to rescale the results back to the original scale.

The resulting average persistence landscapes for degree zero, one and two are in Figure 3, Figure 4, and Figure 5, respectively. The L 1 , L 2 and L ∞ distances between these average landscapes given in Figure 6.

To determine the significance of these distances between average landscapes we performed a permutation test [START_REF] Wasserman | All of statistics[END_REF], which we now explain. First choose a significance level α = 0.05. For every pair i, j ∈ {2, . . . , 10} with i = j, the two corresponding sets of 1000 persistence landscapes were combined in a set of cardinality 2000. Then this set was randomly split into two subsets A 1 , A 2 of cardinality 1000 each. Next the average landscapes λ 1 and λ 2 were computed based on the landscapes in A 1 and A 2 , respectively. Let δ denote the distance between the original average landscapes in dimensions i and j. The distance between λ 1 and λ 2 is compared to δ. This described process is repeated 10, 000 times. The p value equals the proportion of cases in which the distance between λ 1 and λ 2 is greater than δ. For every pair i = j, it never happened that the distance between λ 1 and λ 2 was greater than the corresponding δ. Therefore we can conclude that there is very strong statistical difference between the persistence landscapes in various dimensions.

We believe that the ability to easily perform such calculations will be useful in topological data analysis. . . , 10}. The spheres have been scaled so that the average distance between points is one.

4.2.

Points sampled from [0, 1] d . Given the success of the experiment presented in the Section 4.1 it is natural to ask analogous question for the case of points sampled from a d-dimensional box [0, 1] d . Theorems on the random Vietoris-Rips and Čech complexes can be found in [START_REF] Kahle | Limit theorems for betti numbers of random simplicial complexes[END_REF]. The results presented there indicate that Betti numbers of random complexes "do not live together"; see Figure 1 in [START_REF] Kahle | Limit theorems for betti numbers of random simplicial complexes[END_REF]. Here we focus not on the Betti numbers and their limit distribution, but on a task of dimension detection based on persistent homology. For d ∈ {2, 3, . . . , 10}, to get a single point cloud, we sampled 100 points from [0, 1] d by sampling each coordinate of each of the points independently from an interval [0, 1]. For each dimension d ∈ {2, 3, . . . , 10}, we sampled 1000 point clouds. As in the previous section, to refute the charge that our detected differences are only due to the increase of the average distance between points in higher dimensions, we rescaled the boxes so that the average distance between points in each box is 1. We computed persistent homology in dimension 0 and 1 of the resulting Rips complexes. We obtained 1000 persistence diagrams in dimension 0 and 1 for each dimension d ∈ {2, 3, . . . , 10}.

The permutation test was run for 10, 000 permutations for the obtained persistence intervals separately in dimension 0 and in dimension 1. In all cases, it never happened that the distance between shuffled sets was greater than between the original ones.

This experiment, together with the one presented in the Section 4.1, show that low dimensional persistent homology can be a good tool for estimating the intrinsic dimension of data sets. Knowing the intrinsic dimension of the dataset is important, since when working with point clouds embedded in high dimensional spaces, it is typically assumed that the intrinsic dimension of the data is low and under this assumption one can overcome the curse of dimensionality. The presented techniques may eventually allow one to explicitly estimate the intrinsic dimension of the considered dataset. This was repeated 5 times. For each N the distance between all pairs of persistence diagrams was computed by using the Bottleneck, 1-and 2-Wasserstein metrics and the L ∞ , L 1 and L 2 persistence landscape metric. We used [START_REF] Morozov | The Dionysus software project[END_REF] for the Bottleneck and Wasserstein distance computations. These implementations are known to use sub-optimal/slow algorithms. We are aware of efforts to implement the faster algorithms of [START_REF] Efrat | Geometry helps in bottleneck matching and related problems[END_REF][START_REF] Pankaj | Vertical decomposition of shallow levels in 3-dimensional arrangements and its applications[END_REF], but we are unaware of any that are currently available. We consider the average computation time for each N . The comparison of these times is presented in Figure 7.

4.4.

Computing distance matrices on a random set of persistence intervals. In this experiment, for each N ∈ {100, 200, . . . , 1000}, we generated 1000 random collections of N birthdeath pairs, where the pairs (b, d) chosen independently from the uniform distribution on {0 ≤ b ≤ d ≤ 1} as in Section 4.3. Our aim was to compute the persistence landscape distance matrix for each of these and to see how the computation time scales with the number of birth-death pairs. The computation times are illustrated in the Table 8.

Specification of the implementation

An implementation of the presented procedures is available at the web page http://hans. math.upenn.edu/ ~dlotko/persistenceLandscape.html.

The library is available as a programming tool and can be easily maintained by users who know C++. However our aim is also to provide easy to use tools for users who are not familiar with programming. Therefore a number of programs that use the library have been created for the user's convenience. Also in the future, when there is such demand, we plan to add new programs to this collection. In this section, we describe the programs currently available and illustrate how to use them on a toy example.

These programs have been compiled for Linux, Windows and OS X operating systems, and are also available in the package. All the results described in this Section were obtained by using these programs. More generally, such an input file may encode a linear combination of persistence landscapes. Finally it may be a file containing a list of names of files containing either birth-death pairs or linear combinations of persistence landscapes. Due to its generality, this is a typical input for programs presented in this library. The output files consist of persistence landscapes and linear combinations of persistence landscapes as described above. 5.2. Toy example. Our toy test example for these programs is a set of files consisting of birthdeath pairs, which we now describe. For n ∈ {1, 2, . . . , 5}, let A n denote the union of the circles of radius one centered at (0, 0), (2, 0), . . . , (2n, 0). From A n we sampled 50n points independently using the uniform measure. To each of these points we added a uniform error sampled from [-0.15, 0.15] 2 , see Figure 9 for example. From this set of points, we calculated the persistent homology in degrees zero and one of the corresponding Vietoris-Rips complex using Perseus [START_REF] Mischaikow | Morse theory for filtrations and efficient computation of persistent homology[END_REF]. This was repeated 11 times for each n. The results are encoded in 5 directories (named 1circle,. . . ,5circle) each of which contains 2 × 11 = 22 files that each encode a persistence diagram of either degree zero or one. For the further use in the programs the files containing homology in a particular degree are listed in a files 1circle_dim1.txt,...,5circle_dim1.txt . For example, the file 1circle_dim1.txt lists the following files, each of which contains a persistence diagram. The input parameter is a file containing a list of files containing either persistence diagrams or persistence landscapes from which an average landscape will be calculated. For example, one can call the program computeAverages 1_circle_dim1.txt. As a result, a file containing the average landscape is produced. The same procedure is tested for the remaining files and as a result we obtain ten files with average landscapes.

Plots of landscape(s).

In this section, we demonstrate how to generate output that can be used to plot landscapes. The presented software does not have a built-in graphical engine, so it instead creates gnuplot-readable files. To plot these, the user should install gnuplot [START_REF] Williams | Gnuplot 4.5: an interactive plotting program[END_REF]. The program PlotOfLandscape takes as its first parameter the file name of a file containing either a persistence diagram, a persistence landscape, or a linear combination of persistence landscapes. The remaining two parameters give the range of landscapes which should be plotted. For a range a, b, where a, b ∈ N and a < b the functions λ a , λ a+1 , . . . , λ b-1 will be plotted.

This program was used for each of the average landscapes computed in the Section 5.3 to obtain the plots in Figure 10.

The program PlotsOfLandscapes has as its first parameter the name of a file listing files containing persistence landscapes, with the remaining parameters the same as PlotOfLandscape. It was used to produce Figures 3,4 (1) the name of a file containing names of files containing persistence diagrams, persistence landscapes, or linear combinations of persistence landscapes; and (2) a real p ≥ 1, or p = -1 indicating which norm we want, where p = -1 represents the supremum norm.

The output of the program is a list of the p-norms of the input landscapes.

Here we compute the p-norms of the average landscapes in degree one from Section 5.3. (1) A file with names of files containing persistence diagrams or persistence landscapes or linear combinations persistence landscapes. (2) An integer p. If p ≥ 1, then the L p distance between landscapes will be computed. If p = -1, the L ∞ distance will be computed.

The output is a distance matrix in a text file. So except for the one outlier, the results are exactly as one would expect. 5.9. A classifier based on all dimensions. We have another classifier, ClassifierBasedOnAllDimensions, similar to the one presented in Section 5.8 that uses persistence data from more than one degree.

Exact versus grid based computations.

In this section, we discuss the pros and cons of exact and grid-based computations of persistence landscapes. We do not provide comparison times, since it is not clear what dataset is the right benchmark for such comparisons. For a persistence diagram with a great number of points concentrated in a small region, then a grid-based implementation should be faster. On the other hand, for a persistence diagram in which points appears at a wide range of birth and death parameters, an exact implementation will be superior. This is because representing such a landscape with a reasonable grid spacing (i.e. reasonable accuracy) will require a very large grid. We encourage the users of the Persistence Landscape Toolbox to experiment with both strategies and to pick the better one for the data at hand. To change a software mode from exact computations to a grid-based computations, please modify the self-explanatory file configure in the main folder of the library.

We now discuss error bounds for grid-based estimates of persistence landscapes. The persistence landscape and average persistence landscape are piecewise-linear functions with slopes between -1 and 1. Therefore for each point (x, y) in the plot of a landscape, the landscape will lie between the lines through (x, y) having slope ±1. Now consider two points in the landscape (x 0 , y 0 ) and (x 1 , y 1 ) such that x 0 and x 1 are consecutive points in the grid. Then the persistence landscape between those two points will lie in the intersection of the cones for (x 0 , y 0 ) and (x 1 , y 1 ). These intersecting cones are illustrated in Figure 11. We assume that the grid-based estimate of the persistence landscape is a piecewise linear function through the points obtained by evaluating the persistence landscape on the grid. From the grid and the corresponding values of the persistence landscape one can bound the error. If the grid has spacing δ, then the L ∞ error is bounded by δ 2 . The L 1 error is bounded by half the sum of the areas of the rectangles described above. If the grid has size m and there are K nonzero persistence landscape functions, then this is bounded by K(m -1) δ 2 4 .

  persistence landscape, a list of lists of critical points (x, y). Sort A first according to increasing b and second according to decreasing d; k ← 1; while A = ∅ do Initialize L k ; Pop the first term (b, d) from A; Let p point to the next term; Add (-∞, 0), (b, 0), ( b+d 2 , d-b 2 ) to L k ; while L k . last = (0, ∞) do if d maximal among remaining terms in A starting at p then Add (d, 0), (∞, 0) to L k ; else Let (b , d ) be the first of the terms starting at p with d > d; Pop (b , d ) from A; Let p point to the next term;

Figure 1 .

 1 Figure 1. Algorithm 1 is used to construct the persistence landscape corresponding to the birth-death pairs {(1, 5), (2, 8), (3, 4), (5, 9), (6, 7)}. (a) The functions (3) corresponding to the birth-death pairs and their corresponding critical points. (bd) Steps through the while loop to construct L 1 . (e) The graph of λ 1 . (f) The graph of the functions corresponding to the remaining pairs in the list A. (g-i) The second iteration of the while loop constructs L 2 . (j) The graph of λ 2 . (k) The graphs of the functions corresponding to the remaining pairs on the list A. (l) The graph of λ 3 .

Figure 2 .

 2 Figure2. The worst case scenario for the complexity of the persistence landscapes occurs when all the pairs of persistence intervals have nonempty intersection. Suppose n persistence intervals having nonempty intersection are given. Then there are n nonzero persistence landscapes and the landscape λ k has 2n + 3 -2k critical points. That gives n 2 + 2n = O(n 2 ) critical points in total.

dc

  |ax + b| p dx which can be written as one or two integrals of the form d c (ax + b) p dx = (ax+b) p+1 a(p+1) | d c . This calculation also has time complexity O(P ).

Figure 3 .

 3 Figure3. Average persistence landscapes in degree 0 of points sampled from S d for d ∈ {2, . . . , 10}. The spheres have been scaled so that the average distance between points is one.

Figure 4 .

 4 Figure 4. Average persistence landscapes in degree 1 of points sampled from S d for d ∈ {2, . . . , 10}. The spheres have been scaled so that the average distance between points is one.

Figure 5 . 3 .

 53 Figure 5. Average persistence landscapes in degree 2 of points sampled from S d for d ∈ {2, . . . , 10}. The spheres have been scaled so that the average distance between points is one.

Figure 6 .Figure 7 .

 67 Figure 6. Color plots of L 1 (top), L 2 (middle) and L ∞ (bottom) distance matrices for average persistence landscapes for points sampled from spheres in dimensions {2, . . . , 10} normalized so that average distance between points is 1, for homological degree 0 (left), 1 (middle) and 2 (right).

Figure 8 .

 8 Figure 8. Computation times for the distance matrix computations described in the Section 4.4. Left: the average time to calculated the persistence landscape from a given number of birth-death pairs. Right: the time to calculate the distance matrix of 1000 such persistence landscapes.

Figure 9 .

 9 Figure 9. Sample point cloud from 5-circle directory.

and 5. 5 . 5 .

 55 Norms of landscapes. The program normsOfLandscapes computes norms of persistence landscapes and linear combinations of persistence landscapes. The parameters of the program are:

Figure 10 .

 10 Figure 10. Average persistence landscapes in degree 1 from the toy example in Section 5.2.

1 -

 1 circle 846.959 126.715 27.3333 2-circle 1610.31 174.102 29.0833 3-circle 2462.29 214.864 29.6667 4-circle 3421.15 256.394 31.3333 5-circle 4080.41 278.459 30.5833 5.6. Distance matrix. In this section, we illustrate the usage of the program DistanceMatrix. The parameters of the program are:

Figure 11 .

 11 Figure 11. Exact versus grid-based estimate of the persistence landscape. The underlying persistence landscape is given by the black lines. The blue dots in the x-axis indicate the grid points. The corresponding red dots in the plot of the persistence landscape give the values at the grid points. The dashed line is the estimated persistence landscape. The blue regions indicate the intersection of the cones in which the persistence landscape is guaranteed to be located.

  Average Persistence Landscapes. Let us describe ComputeAverage.

	1circle/points_1circle_example_0_persistence_1.txt
	1circle/points_1circle_example_1_persistence_1.txt
	1circle/points_1circle_example_2_persistence_1.txt
	1circle/points_1circle_example_3_persistence_1.txt
	1circle/points_1circle_example_4_persistence_1.txt
	1circle/points_1circle_example_5_persistence_1.txt
	1circle/points_1circle_example_6_persistence_1.txt
	1circle/points_1circle_example_7_persistence_1.txt
	1circle/points_1circle_example_8_persistence_1.txt
	1circle/points_1circle_example_9_persistence_1.txt
	1circle/points_1circle_example_10_persistence_1.txt
	5.3.

  For example, calling the program DistanceMatrix with the parameters 1_circle_files.txt 2 produces the following. Permutation test. The permutation test was explained in the Section 4.1. Our program PermutationTest performs a permutation test for each pair in a list of lists of persistence diagrams or persistence landscapes, and outputs a matrix of p-values. Its parameters are (5,28.9991) (4,91.9347) (3,154.484) (2,196.559) (1,234.49) (5,38.7583) (4,117.984) (3,172.077) (2,212.901) (1,249.565) (5,44.8902) (4,86.6282) (3,131.766) (2,180.086) (1,220.9) (5,24.9582) (4,105.167) (3,162.33) (2,203.913) (1,238.713)

	0 106.729 161.55 207.872 239.638
	106.729 0 104.46 160.555 199.336
	161.55 104.46 0 105.577 152.975
	207.872 160.555 105.577 0 97.8873
	239.638 199.336 152.975 97.8873 0
	5.7.

Acknowledgments

The authors would like to thank Tamal K. Day, Brittany T. Fasy, Michael Kerber, Miro Kramar, Donald Sheehy and the anonymous referees for their valuable suggestions.

P.B is supported by AFOSR grant FA9550-13-1-0115. P.D is supported by the Advanced Grant of the European Research Council GUDHI (Geometric Understanding in Higher Dimensions), DARPA grant FA9550-12-1-0416 and AFOSR grant FA9550-14-1-0012.

(1) a positive integer M being the number of files with the input persistence intervals or landscapes;

(2) M names of files with the input persistence diagrams or persistence landscapes. Each file is supposed to contain data from the same class, two different files are supposed to contain data from two different classes, (3) Positive integer N indicating the number of tries in the permutation test, (4) An real number p ≥ 1 indicating which distance is to be used in the procedure (or p = -1 for the L ∞ distance). Please note that this procedure usually takes a lot of time. Therefore, after each step a message is outputted to the screen so the user can verify that the program is progressing.

As a sample test we have used all the files with persistence intervals for the N circle for N ∈ {1, . . . , 5}. In this case all p-values are 0 which is what one should expect given how different the persistence landscapes of the different classes are.

5.8.

A classifier based on a single dimension. In this section, present a simple implementation of a nearest-neighbor classifier based on persistence landscapes. This is just one possible example of topological statistics in classification. There are many other ways to perform this task using the software described here. For example, one may apply a support vector machine (SVM) [START_REF] Cortes | Support-vector networks[END_REF] to our distance matrix. Our classifier is implemented in the program ClassifierBasedOnSingleDimension.

The classifier proposed in this paper is based on the following idea. Suppose we are given a training set consisting of N sets of persistence diagrams or persistence landscapes T 1 , . . . , T N . To start, we calculate the average landscape Av 1 , . . . , Av N for each of these classes. Then a sequence of M landscapes l 1 , . . . , l M are given to classification. We give two options:

(1) for l j the program can return an index i ∈ {1, . . . , N } such that the L p distance, for a chosen p, between l j and Av i is the smallest one; or (2) for l j the program can return a vector of pairs (distance, number of class) sorted according to the first coordinate. The usage of the program ClassifierBasedOnSingleDimension is determined by the first parameter which is one of -construct, -classify, or -both.

If the first parameter of the program is -construct, the program will construct the average landscape of each class and write it to a .lan file in the same folder where the program is located. In this case the parameters of the program are:

(1) An integer N indicating the number of classes in the training set, (2) N names of files with each listing files for one of the classes in the training set. In this case the program will write N .lan files with the average landscapes of the N input classes. Those average landscapes can be later used by the program to perform classification when the -classify parameter is used.

If the first parameter of the program is -classify, it is assumed that the average landscapes have already been created (by using -construct option). In this case the parameters of the program are:

(1) A positive integer N indicating how many classes there are in the considered data, (2) The name of a file listing with names of files which are to be classified, (3) A real number p ≥ 1 indicating which norm is to be used in classification (with p = -1 for the supremum norm), (4) A parameter q valued 0 or 1. If q = 0, then for each input persistence diagram or persistence landscape the best matching cluster will be computed. If q = 1, the distances to all averages will be calculated. If the first parameter of the program is -both, then the program both computes the averages of the training set and classifies the test set in one run. In this case the parameters of the program are:

(1) A positive integer N indicating how many classes there are in the considered data, (2) N names of files with each listing files for one of the classes in the training set, (3) The name of file listing the names of the files which are to be classified, (4) A real number p ≥ 1 indicating which norm is to be used in classification, with p = -1 indicating the supremum norm, (5) A parameter q valued 0 or 1 as above.

In the latter two cases the result of the classification is written to the output file classification.txt.

In order to test the classifier we have used the persistence diagrams in degree 1 calculated in Section 5.2. For each class, half of the persistence diagrams were used as a training set and the other half were later classified.

The file listing the files to classify contains the names of all the files which not used in the training set. The files come from 1circle, . . . , 5circle, in order. The results of the classification which outputs only the best matching element are presented below (the formatting has been adjusted to make the interpretation easier). Clearly the classification works very well except from the second case in 4circle where we get the wrong result. To understand the reason for the mismatch we went back to the data and it turned out that one of the intervals corresponding to a circle lived for relatively short time.

When we ask for all the matchings sorted from the best to the worst, we obtain the following.

Even though the slope of an average landscape may be any real number between -1 and 1, in our experience, the slope is often far from the extreme possible values. Therefore the error bounds above may be large overestimates.

Configuration of the library.

The Persistence Landscape Toolbox library is a collection of C++ programs. Each file with the .cpp extension contains a ready-to use program. For a description of its functionality and required parameters, please run the program without any parameters. There is a configure file which needs to be present in the folder in which a program is run. This file contains basic configuration parameters for the library. In this file one can set up the constant which represents infinity. The library assumes that the input file consists of a collection of numbers as described in the Section 5.1. To encode infinite intervals, the infinities have to be changed to a "magic number" which is defined in the configure file. There one can find various options for what can be done with the infinite intervals. This file also allows one to switch the library from the default exact mode to a gridbased mode. If a grid-based mode is used, the user should set up the parameters of the grid in the configuration file.

Conclusion

To conclude, we have provided asymptotically optimal algorithms for computing persistence landscapes, averaging them and calculating distances between (average) persistence landscapes. We have implemented these algorithms and demonstrated how they may be used for hypothesis testing and classification. We hope that they will of use to the (topological) data scientist.

Furthermore, this is not intended to be the end of the story. In future work, we aim to improve these algorithms and related software and increase their usefulness.

Appendix A. Persistence landscapes for barcodes with infinite intervals

We can extend the definition of the persistence landscape in Section 1 to birth-death pairs (b, d) with -∞ ≤ b < d ≤ ∞ using the following definitions.

Algorithm 4 is a generalization of Algorithm 1 that constructs the persistence landscape from a barcode that may contain infinite intervals.

Algorithm 4: Compute the persistence landscape.

Input: A = {(bi, di)} n i=1 -a list of birth-death pairs, -∞ ≤ bi < di ≤ ∞. Output: {L k } -the persistence landscape, a list of lists of critical points (x, y). Sort A first according to increasing b and second according to decreasing d


Romeo Carabelli
email: carabelli@univ-tours.fr

Heritage and sustainable development; bound-partners in a believable future 1

Keywords: Maghreb, Casablanca, Tunis, Algiers, sustainable development, modern heritage

Two global notions are improving their importance: heritage and sustainable developmentobviously involving safeguard and climate changes. New transformations in urban areas are showing the great tension existing among heritage and development needs. Heritage and Sustainable development still need, public interventions with the purpose of implement a virtuous and self standing development strategy. The great Maghreb's towns will be presented as study case: challenges, chances and realities.

It's largely clear that the main issues to manage heritage and sustainability at urban scale are not any more technical ones. The issue is a global policy able to combine together the countless partial strategies and interests existing in human settlements.

Generally speaking, both heritage and sustainability are requesting public strategies and policies; therefore a transformation of general urban actors in advised ones is needed.

Paying attention to the general speeches about town planning, urban transformation and landscape modification, we can frequently observe the use, and often the abuse, of two global notions: heritage and sustainable developmentobviously involving safeguard and climate changes.

Those two notions are global ones since they are inclusive and each one of them is ready to contain the whole ecumene inside its own universal space. Of course I'm overstating the reality, but not so much. Actually these two "expanding" notions will be obliged to face each other and they are called to manage a common way of life.

Will they find a pacific co-existence or will they suffer of an obliged one?Will the future of our globe assist to a competition between sustainability(s) and heritage(s) or they will collaborate among them?

We are interested in following the development of the relations among those different ways to face the contemporary changes mainly because the two notions look like based over a common acting area. They are overlapping at least on two fields we are really interested in; both are summon up intangible assets in order to justify a way of working on tangible ones. They are also referring to an extended time, largely wider then a direct personal experience; we like to imagine a special site on the timeline that's almost "past of past" and "future of future". In that special site, the two notions are collaboratingor not.

Inflecting sustainable development into urban landscapes transformation, we can learn about ecoquarters, urban project involved in global sustainability, passive and active houses or, more generically, sustainable habitat. Those actions are really interesting and central in a new way of planning and conceptualizing the nowadays spaces for human life. Heritage and Development Cooperation: 16th, 17th and 18th June, 2010. Seville (Spain)., Sevilla, Comité Científico del IV Congreso de Patrimonio Cultural y Cooperación al Desarrollo, 2010, pp. 85-89. At the same time, we know that more of 50% of world inhabitants are already living in urban areas; not only in developed country but also in developing ones. It means that the interest on new awareness about climate change and energy issues is crucial but … but it need to be supported by a special issue on urban built environment, the real estate environment for the most of human beings.

All new sustainable building, all new sustainable quarter is a laudable one but it doesn't solve the problem if we don't manage the (old) existing habitat. People going to live in a "green", "low impact" and "sustainable" environment; leave their former environment, which is not changing. Therefore, the new sustainability is forgetting the old estate, which will be used by other people, improving the global environmental impact.

Should we destroy our built environment if is not sustainable? Maybe yes, maybe no but in some parts of our towns, the question is a displaced one at least because the heritage areas are areas we decided to preserve. In some urban places, the two notions must manage a way to live together because heritage is there and climate is changing. How to manage it? Will be a cultural synthesis between the two notions possible? Architectural and urban patrimonial inheritance are components of ours life, therefore they have to combine theirs transformations integrating living societies; they are a part of the contemporary existence, consequently they are involved in the resiliency action towards the new sustainable needs.

To proceed with field observation, we have chosen to work in the specific area of Mediterranean towns because of a global reason: they are both subjected to a "patrimonialisation" process and immediately connected with some climate changes, and more precisely with the climate heating and the water rarefaction.

The historical studies and the practical observation in different regions of the Mediterranean highlight the existence of a "new" architectural and urban modern heritage dating from the 19th and the 20th centuries, largely unknown and underestimated. This heritage results from the mutual interaction of different cultures and the merging of local and exogenous know-how that occurred among the European colonies, protectorates and mandates and in the late Ottoman Empire but it was not only an outcome of the "western and colonial" rule, as it is frequently considered.In many cases indeed, the southern city fabrics were also invested by modernization processes, integrating new technologies and spatial configurations in the "traditional" urban structure.

New types of buildings and urban spaceswere created, which are appropriated and used by the present urban societies, with a specific richness and meaning. Their conservation requires the understanding of multiple and often conflicting models and meanings, but also theraising of the awareness and the dissemination of their mutual contribution and interpenetration in the different Mediterranean cultures.

Architectural and urban features of the last two centuries are nowadays spatial markers, which are shaping our towns. This material production is associated to the industrial revolution, to the general modernisationof all the Mediterranean countries, to some national "westernalization" strategies and, of course, to colonialism.

At the same time, those quarters are at the core of the nowadays urban areas and they need to be transformed in order to participate to the contemporary challenge of economical, social and sustainable development. Therefore, all countries are keeping into their agenda new changes to update their towns.

Those recent quarters are indicated as "new towns" (villesnouvelles), or "European towns", or "colonial towns" to differentiate them from the ancient walled cities, the medinas as they are nowadays called in Maghreb. All those definitions are correct and incorrect at the same time, as all stereotypes they are an absolute simplification of the complex meaning of a human and social realities. There were new parts of the town also inside the ancient quarters; while outside the walls there were also local/native inhabitants, and not only European or colonial ones.

These modern areas are interesting for our activities also because of theirs typologies; they are modern and active quarters, free of old fashioned patrimonial identities. The potential value of the last two centuries architectural and urban heritage needs to be enhanced and requires a better valorisation to play a proactive role in the development strategies, at the national and local level.

These 19 th and 20 th century urban realities often reflect specific and relevant social and economical configurations, in spite of the lack of appropriate conservation and development initiatives respecting their tangible and intangible heritage.

We estimate that the mutual heritage of the last two centuries across Mediterranean basin represents a common field of action that needs a global approach including a significant interest to sustainable methods of transformation.Due to its recentand often imported and imposedorigin, this heritage, despite its massive presence in the historic fabric and townscape, is rather neglected and suffers from a lack of interest; clearly reflected in the inadequacy of valorisation strategies and sensitization campaigns.

New transformations in those areas are showing the great tension existing among heritage and development needs. In this context, the only chance we have to propose to next generations a safeguarded heritage included into a sustainable development strategy is to enhance knowledge and awareness about heritage and sustainability and, at the same time, to improve public leaded strategies.

Heritage needed, and still needs, public interventions with the purpose of implement safeguard. Sustainable development too needs public strategies to take off towards a virtuous and self standing implementing strategy. Heritage and sustainability are collective and global interests; the private sector is frequently unable to manage global and collective issues, it has not the intellectual and operational instruments to deal with this kind of stakes.

The goal for a winning strategy should be the transformation of the diffuse heritage into contemporary active and participation instruments of development. The great Maghreb's towns will be presented as study case: challenges, chances and realities.

We realized that the ground level is the richest place for active interfaces among the town's formal and informal actors and the main place to deal with. That's the reason why we are mainly focusing on this ground level, the "0" one, the one of the public road.

The three biggest towns in Maghreb area -Algiers, Casablanca and Tunisare under transformation through urban projects that are producing spin offs in both areas of sustainable development and heritage safeguard. Obviouslythey are producing different kind of heritage areas and, at the same time, they are introducing indispensable urban innovations, perhaps sustainable ones.

Tunis downtown is an old laboratory case; in order to avoid all unsuitable references to the colonial period the area is called "hyper centre". It's a space which includes nowadays the whole ancient walled city and the more representative area built during the last two centuries. With or without consciousness, the "heritage core" is becoming anunicum, from the "lake" to the Kasbah's hill.

We point out that this large area was already conceptualized in the late 1960s by the overall project proposed by the Olivier-Clément Cacoub atelier during the first period of Bourquiba presidency. That project, forecasting a global transformation of the area following the International style concepts, conceived wide and irresponsible demolitions in the old town and, luckily, it wasn't achieved.

Several works of heritage protection were fulfilled in Tunis and in the late 1990s the rehabilitation of the Bourguiba Avenue was the chance to (re)establish the conceptual cohesion of the spaces between the lake and the hill. The process was really controversial and the result is really different from the initially forecast configuration.

In Algiers, after some decades of slow motion following political and terrorism issues, the 19 th and 20 th century town is dealing with a wide programme of heavy infrastructures of public transport. The whole "recent historical centre" (the historical centre without the well-known Kasbah) will be encompassed by rail transports, tramways and metro. This will transform completely the approach to the recent heritage centre and to the private transport system.

Men are at work everywhere in Casablanca's nuke. After the pharaonic project of the great mosque, the coast is "on work": a new look for the Corniche, the erection of a high class Marina and the reconstruction of CasaPort train station area. Immediately inside those areas, a tramway service is foreseen, even if the project had languishing during the last two years, at the opposite of the tramway in the capital, Rabat, where the tramway line has quickly improved.

It's largely clear that the main issues to manage heritage and sustainability at urban scale are not any more technical ones; nowadays technology and experts exist on site, locals or internationals, while once the lack of expert and technology was a real handicap for both heritage safeguard and sustainable development. The issue is a global one able to combine together the countless partial strategies and interests existing in a town and formulating global policies for an adequate development.

Generally speaking, both heritage and sustainability are requesting public strategies and policies. Only with a global urban approachit is possible todesign solutions for those ambitions;partial actionsfocusing only heritage or only sustainability can't be efficient ones. Actually, the first step to improve a positive operative merge between safeguard and sustainability is a friendly analysis system able to offer an appropriate environment to "natural" actors. Fields operations need to transform ordinary landscape and urban actors in advised ones.

BIBLIOGRAPHICAL REFERENCES

This text is a partial formalisation of two ongoing experiences in the frame of two projects that the author is coordinating:

1. Mutual Heritage: from historical integration to contemporary active participation, a project on the recent architectural and urban heritage in the Mediterranean area funded by the European Union within the Euromed Heritage 4 programme.

1

  CARABELLI, ROMEO, « Heritage and sustainable development; bound-partners in a believable future », inCONGRESO INTERNACIONAL DE PATRIMONIO CULTURAL Y COOPERACIÓN AL DESARROLLO (ed.), Libro de actas del IV Congreso Internacional de Patrimonio Cultural y Cooperación al Desarrollo 16, 17 y 18 de junio de 2010, Sevilla (España) = Prints of the papers of the IV International Congress on Cultural

2.

  R+0 !Développement durable et conception des espaces publics des centres modernes des villes méditerranéennes a project on sustainabledevelopment in recentheritagespacesfunded by the PACA region, the PUCA and the Caisse des dépots et des consignations within the Projet d'Institut du Bâtiment Méditerranéen (IBMed)