Persistent Cohomology and Circular Coordinates
Abstract
Nonlinear dimensionality reduction (NLDR) algorithms such as Isomap, LLE and Laplacian Eigenmaps address the problem of representing highdimensional nonlinear data in terms of lowdimensional coordinates which represent the intrinsic structure of the data. This paradigm incorporates the assumption that realvalued coordinates provide a rich enough class of functions to represent the data faithfully and efficiently. On the other hand, there are simple structures which challenge this assumption: the circle, for example, is onedimensional but its faithful representation requires two real coordinates. In this work, we present a strategy for constructing circlevalued functions on a statistical data set. We develop a machinery of persistent cohomology to identify candidates for significant circlestructures in the data, and we use harmonic smoothing and integration to obtain the circlevalued coordinate functions themselves. We suggest that this enriched class of coordinate functions permits a precise NLDR analysis of a broader range of realistic data sets.
Persistent Cohomology and Circular Coordinates
Vin de Silva 
Department of Mathematics 
Pomona College 
Claremont, California 
vin.desilva@pomona.edu 
Mikael VejdemoJohansson 
Department of Mathematics 
Stanford University 
Stanford, California 
mik@math.stanford.edu 
\@float
copyrightbox[b]
\end@floatCategories and Subject Descriptors G.3 [Probability and statistics]: multivariate statistics; I.5.1 [Pattern recognition]: models—geometric

algorithms, theory

dimensionality reduction, computational topology, persistent homology, persistent cohomology
Nonlinear dimensionality reduction (nldr) algorithms address the following problem: given a highdimensional collection of data points , find a lowdimensional embedding (for some ) which faithfully preserves the ‘intrinsic’ structure of the data. For instance, if the data have been obtained by sampling from some unknown manifold — perhaps the parameter space of some physical system — then might correspond to an dimensional coordinate system on . If is completely and nonredundantly parametrized by these coordinates, then the nldr is regarded as having succeeded completely.
Principal components analysis, or linear regression, is the simplest form of dimensionality reduction; the embedding function is taken to be a linear projection. This is closely related to (and sometimes identifed with) classical multidimensional scaling [?].
When there are no satisfactory linear projections, it becomes necessary to use nldr. Prominent algorithms for nldr include Locally Linear Embedding [?], Isomap [?], Laplacian Eigenmaps [?], Hessian Eigenmaps [?], and many more.
These techniques share an implicit assumption that the unknown manifold is welldescribed by a finite set of coordinate functions . Explicitly, some of the correctness theorems in these studies depend on the hypothesis that has the topological structure of a convex domain in some . This hypothesis guarantees that good coordinates exist, and shifts the burden of proof onto showing that the algorithm recovers these coordinates.
In this paper we ask what happens when this assumption fails. The simplest space which challenges the assumption is the circle, which is onedimensional but requires two real coordinates for a faithful embedding. Other simple examples include the annulus, the torus, the figure eight, the 2sphere, the last three of which present topological obstructions to being embedded in the Euclidean space of their natural dimension. We propose that an appropriate response to the problem is to enlarge the class of coordinate functions to include circlevalued coordinates . In a physical setting, circular coordinates occur naturally as angular and phase variables. Spaces like the annulus and the torus are well described by a combination of real and circular coordinates. (The 2sphere is not so lucky, and must await its day.)
The goal of this paper is to describe a natural procedure for constructing circular coordinates on a nonlinear data set using techniques from classical algebraic topology and its 21stcentury grandchild, persistent topology. We direct the reader to [?] as a general reference for algebraic topology, and to [?] for a streamlined account of persistent homology.
There have been other attempts to address the problem of finding good coordinate representations of simple nonEuclidean data spaces. One approach [?] is to use modified versions of multidimensional scaling specifically devised to find the best embedding of a data set into the cylinder, the sphere and so on. The target space has to be chosen in advance. Another class of approaches [?, ?] involves cutting the data manifold along arcs and curves until it has trivial topology. The resulting configuration can then be embedded in Euclidean space in the usual way. In our approach, the number of circular coordinates is not fixed in advance, but is determined experimentally after a persistent homology calculation. Moreover, there is no cutting involved; the coordinate functions respect the original topology of the data.
The principle behind our algorithm is the following equation from homotopy theory, valid for topological spaces with the homotopy type of a cell complex (which covers everything we normally encounter):
(1) The lefthand side denotes the set of equivalence classes of continuous maps from to the circle ; two maps are equivalent if they are homotopic (meaning that one map can be deformed continuously into the other); the righthand side denotes the 1dimensional cohomology of , taken with integer coefficients. In other language: is the classifying space for , or equivalently is the Eilenberg–MacLane space . See section 4.3 of [?].
If is a contractible space (such as a convex subset of ), then and Equation (Persistent Cohomology and Circular Coordinates) tells us not to bother looking for circular functions: all such functions are homotopic to a constant function. On the other hand, if has nontrivial topology then there may well exist a nonzero cohomology class ; we can then build a continuous function which in some sense reveals .
Our strategy divides into the following steps.

Represent the given discrete data set as a simplicial complex or filtered simplicial complex.

Use persistent cohomology to identify a ‘significant’ cohomology class in the data. For technical reasons, we carry this out with coefficients in the field of integers modulo , for some prime . This gives us .

Lift to a cohomology class with integer coefficients: .

Smoothing: replace the integer cocycle by a harmonic cocycle in the same cohomology class: .

Integrate the harmonic cocycle to a circlevalued function .
The paper is organized as follows. In Section Persistent Cohomology and Circular Coordinates, we derive what we need of equation (Persistent Cohomology and Circular Coordinates). Steps (1–5) of the algorithm are addressed in Sections Persistent Cohomology and Circular Coordinates–Persistent Cohomology and Circular Coordinates, respectively. In Section Persistent Cohomology and Circular Coordinates we report some experimental results.
Let be a finite simplicial complex. Let denote the sets of vertices, edges and triangles of , respectively. We suppose that the vertices are totally ordered (in an arbitrary way). If then the edge between vertices is always written and not . Similarly, if then the triangle with vertices is always written .
Cohomology can be defined as follows. Let be a commutative ring (for example ). We define 0cochains, 1cochains, and 2cochains as follows:
These are modules over . We now define coboundary maps and .
Let . If we say that is a cocycle. If admits a solution we say that is a coboundary. The solution , if it exists, can be thought of as the discrete integral of . It is unique up to adding constants on each connected component of .
It is easily verified that for any . Thus, coboundaries are always cocycles, or equivalently . We can measure the difference between coboundaries and cocycles by defining the 1cohomology of to be the quotient module
We say that two cocycles are cohomologous if is a coboundary.
We now consider integer coefficients. The following proposition fulfils part of the promise of equation (Persistent Cohomology and Circular Coordinates), by producing circlevalued functions from integer cocycles. It will be helpful to think of as the quotient group .
Proposition 1
Let be a cocycle. Then there exists a continuous function which maps each vertex to 0, and each edge around the entire circle with winding number .

We can define inductively on the vertices, edges, triangles, … of missing. The vertices and edges follow the prescription in the statement of the proposition. To extend to the triangles, it is necessary that the winding number of along the boundary of each triangle is zero. And indeed this is . Since the higher homotopy groups of are all zero ([?], section 4.3), can then be extended to the higher cells of without obstruction.
The construction in Proposition 1 is unsatisfactory in the sense that all vertices are mapped to the same point. All variation in the circle parameter takes place in the interior of the edges (and higher cells). This is rather unsmooth. For more leeway, we consider real coefficients.
Proposition 2
Let be a cocycle. Suppose we can find and such that . Then there exists a continuous function which maps each edge linearly to an interval of length , measured with sign.
In other words, we can construct a circlevalued function out of any real cocycle whose cohomology class lies in the image of the natural homomorphism .

Define on the vertices of missing by setting to be mod . For each edge , we have
which is congruent to mod , since is an integer.
It follows that can be taken to map linearly onto an interval of signed length . Since is a cocyle, can be extended to the triangles as before; then to the higher cells.
Proposition 2 suggests the following tactic: from an integer cocycle we construct a cohomologous real cocycle , and then define mod on the vertices of . If we can construct so that the edgelengths are small, then the behaviour of will be apparent from its restriction to the vertices. See Section Persistent Cohomology and Circular Coordinates.
We now begin describing the workflow in detail. The input is a pointcloud data set: in other words, a finite set or more generally a finite metric space. The first step is to convert into a simplicial complex and to identify a stablelooking integer cohomology class. This will occupy the next three subsections.
The first lesson of pointcloud topology [?] is that pointclouds are best represented by 1parameter nested families of simplicial complexes. There are several candidate constructions: the Vietoris–Rips complex has vertex set and includes a simplex whenever all vertices lie pairwise within distance of each other. The witness complex uses a smaller vertex set and includes a simplex when the vertices lie close to other points of , in a certain precise sense (see [?, ?]). In both cases, whenever . Either of these constructions will serve our purposes, but the witness complex has the computational advantage of being considerably smaller.
We determine only up to its 2skeleton, since we are interested in .
Having constructed a 1parameter family , we apply the principle of persistence to identify cocycles that are stable across a large range for . Suppose that are the critical values where the complex gains new cells. The family can be represented as a diagram
of simplicial complexes and inclusion maps. For any coefficient field , the cohomology functor converts this diagram into a diagram of vector spaces and linear maps over ; the arrows are reversed:
According to the theory of persistence [?, ?], such a diagram decomposes as a direct sum of 1dimensional terms indexed by halfopen intervals of the form . Each such term corresponds to a cochain that satisfies the cocycle condition for and becomes a coboundary for . The collection of intervals can be displayed graphically as a persistence diagram, by representing each interval as a point in the Cartesian plane above the main diagonal. We think of long intervals as representing trustworthy (i.e. stable) topological information.
Choice of coefficients. The persistence decomposition theorem applies to diagrams of vector spaces over a field. When we work over the ring of integers , however, the result is known to fail: there need not be an interval decomposition. This is unfortunate, since we require integer cocycles to construct circle maps. To finesse this problem, we pick an arbitrary prime number (such as ) and carry out our persistence calculations over the finite field . The resulting cocyle must then be converted to integer coefficients: we address this in Section Persistent Cohomology and Circular Coordinates.
In principle we can use the ideas in [?] to calculate the persistent cohomology intervals and then select a long interval and a specific . We then let and take to be the cocycle in corresponding to the interval.
Explicitly, persistent cocycles can be calculated in the following way. We thank Dmitriy Morozov for this algorithm. Suppose that the simplices in the filtered complex are totally ordered, and labelled so that arrives at time . For we maintain the following information:

a set of indices associated with ‘live’ cocycles;

a list of cocycles in .
The cocycle involves only and those simplices of the same dimension that appear later in the filtration sequence (thus only with ).
Initially and the list of cocycles is empty.
To update from to , we compute the coboundaries of the cocycles of within the larger complex obtained by including the simplex . In fact, these coboundaries must be multiples of the elementary cocycle defined by , and otherwise. We can write . If all the are zero, then we have one new cocycle: let and define . Otherwise, we must lose a cocycle. Let be the largest index for which . We delete by setting , and we restore the earlier cocycles by setting . In this latter case, we write the persistence interval to the output.
At the end of the process, surviving cocycles are associated with semiinfinite intervals: for .
Remark. The reader may be more familiar with persistence diagrams in homology rather than cohomology. In fact, the universal coefficient theorem [?] implies that the two diagrams are identical. The salient point is that cohomology is the vectorspace dual of homology, when working with field coefficients. That said, we cannot simply use the usual algorithm for persistent homology: we are interested in obtaining explicit cocycles, whereas the classical algorithm [?] returns cycles.
We will establish the correctness of this algorithm in the archival version of this paper. The expert reader may regard this as an exercise in the theory of persistence.
We now have a simplicial complex and a cocycle . The next step is to ‘lift’ by constructing an integer cocycle which reduces to modulo .
To show that this is (almost) always possible, note that the short exact sequence of coefficient rings gives rise to a long exact sequence, called the Bockstein sequence (see Section 3.E of [?]). Here is the relevant section of the sequence:
By exactness, the Bockstein homomorphism induces an isomorphism between the cokernel of and the kernel of , and this kernel is precisely the set of torsion elements of . If there is no torsion, then it follows immediately that the cokernel of the first map is zero. In other words is surjective; any cocycle can be lifted to a cocycle .
If we are unluckily sabotaged by torsion, then we pick another prime and redo the calculation from scratch: it is enough to pick a prime that does not divide the order of the torsion subgroup of , so almost any prime will do.
In practice, we construct by taking the coefficients of in and replacing them with integers in the correct congruence class modulo . The default choice is to choose coefficients close to zero. If then we are done; otherwise it becomes necessary to do some repair work. Certainly modulo , so we can write for some . In the absence of torsion, we can then solve for , and then the required lift is . Fortunately, this has not proved necessary in any of our examples.
Remark. We expect that torsion is extremely rare in ‘real’ data sets, since it is symptomatic of rather subtle topological phenomena. For instance, the simplest examples which exhibit 2torsion are the nonorientable closed surfaces (such as the projective plane and the Klein bottle).
Given an integer cocycle , or indeed a real cocycle , we wish to find the ‘smoothest’ real cocycle cohomologous to . It turns out that what we want is the harmonic cocycle representing the cohomology class .
We define smoothness. Each of the spaces comes with a natural Euclidean metric:
A circlevalued function is ‘smooth’ if its total variation across the edges of is small. The terms capture the variation across individual edges; therefore what we must minimize is .
Proposition 3
Let . There is a unique solution to the leastsquares minimization problem
(2) Moreover, is characterized by the equation , where is the adjoint of with respect to the inner products on .

Note that if then for any we have
which implies that such an must be the unique minimizer. For existence, note that
certainly has a solution missing if . But this is a standard fact in finitedimensional linear algebra: for any real matrix ; this follows from the singular value decomposition, for instance.
Remark. It is customary to construct the Laplacian . The twin equations and immediately imply (and conversely, can be deduced from) the single equation ; in other words is harmonic.
The leastsquares problem in equation (2) can be solved using a standard algorithm such as LSQR [?]. By Proposition 2 we can use the solution parameter to define the circular coordinate on the vertices of . This works because the original cocycle has integer coefficients.
More generally, if is an arbitrary real cocycle such that , it is a straightforward matter to integrate to a circlevalued function on the vertex set . Suppose that is connected (if not, each connected component can be treated separately) and pick a starting vertex and assign . One can use Dijkstra’s algorithm to find shortest paths to each remaining vertex from . When a new vertex enters the structure via an edge , we assign (or if the edge is correctly identified as ). If a vertex is connected to by multiple paths then the different possible values of differ by an integer; this is where we use the hypothesis that is cohomologous to an integer cocyle.
The following experiments were carried out using the Javabased jPlex simplicial complex software [?], with highlevel scripting and numerical analysis in MATLAB. We ran a development version of jPlex to obtain explicit persistent cohomology cocycles. We expect to include the code in the next release of jPlex. We used Paige and Saunders’ implementation of LSQR [?] for the leastsquares problem in the harmonic smoothing step. Timings were determined using MATLAB’s builtin ‘tic’ and ‘toc’ commands, and are included for relative comparison against each other.
We tested our methods on several synthetic data sets with known topology, ranging from the humble circle itself to a genus2 surface (‘double torus’). Most of the examples were embedded in or , with the exception of a sample from a complex projective curve (embedded in ) and a synthetic imagelike data set (embedded in ).
In each case we selected vertices for the filtered simplicial complex: either the whole set, or a smaller welldistributed subset of ‘landmarks’ selected by iterative furthestpoint sampling. We then built a Rips or witness complex, with maximum radius generally chosen to ensure around simplices in the complex.
In most cases, we show the persistence diagram produced by the cocycle computation. The chosen value is marked on the diagonal, with its upperleft quadrant indicated in green lines. The persistent cocycles available at that parameter value are precisely those contained in that quadrant. Each of those cocycles then produces a circular coordinate.
There are various figures associated with each example. Most important are the correlation scatter plots: each scatter plot compares two circular coordinate functions. These may be functions produced by the computation (‘inferred coordinates’) or known parameters. These scatter plots are drawn in the unit square, which is of course really a torus .
When the original data are embedded in or , we also display the circular coordinates directly on the data set, plotting each point in color according to its coordinate value interpreted on the standard huecircle. This works less well in grayscale reproductions, of course.
Finally, in certain cases we plot coordinate values against frequency, as a histogram. This distributional information can sometimes be useful in the absence of other information.
Remark. When the goal is to infer the topology of a data set whose structure is unknown, we do not have any ‘known parameters’ available to us. We can still construct correlation scatter plots between pairs of inferred coordinates, and the distributional histograms for each coordinate individually. We exhort the reader to view the following examples through the lens of the topological inference problem: what structures can be distinguished using scatter plots and histograms (and persistence diagrams) alone?
We begin with the circle itself, and its tautological circlevalued coordinate.
We picked 400 points distributed along the unit circle. We added a uniform random variable from to each coordinate. A Rips complex was constructed with maximal radius 0.5, resulting in 23475 simplices. The computation of cohomology finished in 237 seconds.
Parametrizing at 0.4 yielded a single coordinate function, which very closely reproduces the tautological angle function. Parametrizing at 0.14 yielded several possible cocycles. We selected one of those with low persistence; this produced a parametrization which ‘snags’ around a small gap in the data.
See Figure Persistent Cohomology and Circular Coordinates.
The left panel in each row shows the histogram of coordinate values; the middle panel shows the correlation scatter plot against the known angle function; the right panel displays the coordinate using color. The highpersistence (‘global’) coordinate correlates with the angle function with topological degree 1. Variation in that coordinate is uniformly distributed, as seen in the histogram. In contrast, the lowpersistence (‘local’) coordinate has a spiky distribution.
Another example with circle topology: see Figure Persistent Cohomology and Circular Coordinates. We picked 400 points distributed along the torus knot on a torus with radii 2.0 and 1.0. We jittered them by a uniform random variable from added to each coordinate. We generated a Rips complex up to radius 1.0, acquiring 36936 simplices. We computed persistent cohomology in 70 seconds.
Figure \thefigure: Trefoil torus knot. Persistence diagram (left), correlation scatter plot of inferred coordinate against known parametrization (middle), inferred coordinate in color (right). As expected, the inferred coordinate correlates strongly with the known parameter with topological degree 1. The histogram shows three ‘bulges’ corresponding to the three highdensity regions of the sampled curve, which occur when the curve approaches the central axis of the torus.
For a more elaborate data set with topology, we generated a sequence of 657 rendered images of a colorful cube rotating around one axis. Each image was regarded as a vector in the Euclidean space . From this data we built a witness complex with 50 landmark points and constructed a single circular coordinate. Interpolating the resulting function linearly between the landmarks gave us coordinates for all the points in the family.
See Figure Persistent Cohomology and Circular Coordinates. The frequency distribution is comparatively smooth (by which we mean that there are no large spikes in the histogram), which indicates that the coordinate does not have large static regions. The correlation plot of the inferred coordinate against the original known sequence of the cube images shows a correlation with topological degree 1. We show the progression of the animation on an evenlyspaced sample of representative points around the circle.
See Figure Persistent Cohomology and Circular Coordinates for these two examples.
Conjoined circles: we picked 400 points distributed along circles in the plane with radius 1 and with centres at . The points were then jittered by adding noise to each coordinate taken uniformly randomly from the interval . A Rips complex was constructed with maximal radius 0.5, resulting in 76763 simplices. The cohomology was computed in 378 seconds.
Disjoint circles: 400 points were distributed on circles of radius 1 centered around in the plane. These points were subsequently disturbed by a uniform random variable from . We constructed a Rips complex with maximum radius 0.5, which gave us 45809 simplices. The cohomology computation finished in about 117 seconds.
In both cases, our method detects the two most natural circlevalued functions. The scatter plots appear very similar. In the conjoined case, there is some interference between the two circles, near their meeting point.
See Figure Persistent Cohomology and Circular Coordinates. We picked 400 points at random in the unit square, and then used a standard parametrization to map the points onto a torus with inner and outer radii 1.0 and 3.0. These were subsequently jittered by adding a uniform random variable from to each coordinate. We constructed a Rips complex with maximal radius , resulting in 61522 simplices. The corresponding cohomology was computed in 209 seconds.
Figure \thefigure: Torus in . The two inferred coordinates in this (fairly typical) experimental run recover the original coordinates essentially perfectly: the first inferred coordinate correlates with the meridional coordinate with topological degree , while the second inferred coordinate correlates with the longitudinal coordinate with degree .
When the original coordinates are unavailable, the important figure is the inferredversusinferred scatter plot. In this case the scatter plot is fairly uniformly distributed over the entire coordinate square (i.e. torus). In other words, the two coordinates are decorrelated. This is slightly truer (and more clearly apparent in the scatter plot) for the two original coordinates. Contrast these with the corresponding scatter plots for a pair of circles (conjoined or disjoint).
See Figure Persistent Cohomology and Circular Coordinates. For fun, we repeated the previous experiment with a torus abstractly defined as the zero set of a homogeneous cubic polynomial in three variables, interpreted as a complex projective curve. We picked 400 points at random on , subject to the cubic equation
To interpret these as points in , we used the projectively invariant metric
for all pairs . With this metric we built a Rips complex with maximal radius 0.15. The resulting complex had 44184 simplices, and the cohomology was computed in 56 seconds. We found two dominant coclasses that survived beyond radius 0.15, and we computed our parametrizations at the 0.15 mark.
The resulting correlation plot quite clearly exhibits the decorrelation which is characteristic of the torus.
Figure \thefigure: Elliptic curve. Persistence diagram (left), correlation scatter plot between the two coordinates (right). See Figure Persistent Cohomology and Circular Coordinates. We constructed a genus2 surface by generating 1600 points on a torus with inner and outer radii 1.0 and 3.0; slicing off part of the data set by a plane at distance 3.7 from the axis of the torus, and reflecting the remaining points in that plane. The resulting data set has 3120 points. Out of these, we pick 400 landmark points, and construct a witness complex with maximal radius 0.6. The landmark set yields a covering radius and a complex with 70605 simplices. The computation took 748 seconds active computer time. We identified the four most significant cocycles.
Figure \thefigure: Double torus in . Note that coordinates 1 and 4 are ‘coupled’ in the sense that they are supported over the same subtorus of the double torus. The scatter plot shows that the two coordinates appear to be completely decorrelated except for a large mass concentrated at a single point. This mass corresponds to the other subtorus, on which coordinates 1 and 4 are essentially constant. A similar discussion holds for coordinates 2 and 3.
The uncoupled coordinate pairs (1,2), (1,3), (2,4), (3,4) produce scatter plots reminiscent of two conjoined or disjoint circles.
We are immensely grateful to Dmitriy Morozov: he has given us considerable assistance in implementing the algorithms in this paper. In particular we thank him for the persistent cocycle algorithm.
Thanks also to Jennifer Kloke for sharing her analysis of a visual image data set; this example did not make the present version of this paper. Finally, we thank Gunnar Carlsson, for his support and encouragement as leader of the topological data analysis research group at Stanford; and Robert Ghrist, as leader of the DARPAfunded project Sensor Topology and Minimal Planning (SToMP).
 1 M. Belkin and P. Niyogi. Laplacian Eigenmaps and spectral techniques for embedding and clustering. In T. Diettrich, S. Becker, and Z. Ghahramani, editors, Advances in Neural Information Processing Systems 14, pages 585–591. MIT Press, Cambridge, Massachussetts, 2002.
 2 T. F. Cox and M. A. A. Cox. Multidimensional Scaling. Chapman & Hall, London, 1994.
 3 V. de Silva and G. Carlsson. Topological estimation using witness complexes. In M. Alexa and S. Rusinkiewicz, editors, Eurographics Symposium on PointBased Graphics, ETH, Zürich, Switzerland, 2004.
 4 M. Dixon, N. Jacobs, and R. Pless. Finding minimal parameterizations of cylindrical image manifolds. In CVPRW ’06: Proceedings of the 2006 Conference on Computer Vision and Pattern Recognition Workshop, page 192, Washington, DC, USA, 2006. IEEE Computer Society.
 5 D. L. Donoho and C. Grimes. Hessian Eigenmaps: New locally linear embedding techniques for highdimensional data. Technical Report TR 200308, Department of Statistics, Stanford University, 2003.
 6 H. Edelsbrunner, D. Letscher, and A. Zomorodian. Topological persistence and simplification. Discrete and Computational Geometry, 28:511–533, 2002.
 7 R. Ghrist. Barcodes: the persistent topology of data. Bulletin of the American Mathematical Society, 45(1):61–75, 2008.
 8 L. J. Guibas and S. Y. Oudot. Reconstruction using witness complexes. In Proc. 18th ACMSIAM Sympos. on Discrete Algorithms, pages 1076–1085, 2007.
 9 A. Hatcher. Algebraic Topology. Cambridge University Press, Cambridge, 2002.
 10 J. A. Lee and M. Verleysen. Nonlinear dimensionality reduction of data manifolds with essential loops. Neurocomputing, 67:29–53, 2005.
 11 C. C. Paige and M. A. Saunders. LSQR: Sparse equations and least squares. http://www.stanford.edu/group/SOL/software/lsqr.html.
 12 C. C. Paige and M. A. Saunders. LSQR: An algorithm for sparse linear equations and sparse least squares. ACM Transactions on Mathematical Software, 8(1):43–71, March 1982.
 13 R. Pless and I. Simon. Embedding images in nonflat spaces. In Conference on Imaging Science Systems and Technology, pages 182–188, 2002.
 14 S. Roweis and L. Saul. Nonlinear dimensionality reduction by locally linear embedding. Science, 290:2323–2326, Dec. 2000.
 15 H. Sexton and M. VejdemoJohansson. jPlex simplicial complex library. http://comptop.stanford.edu/programs/jplex/.
 16 J. B. Tenenbaum, V. de Silva, and J. C. Langford. A global geometric framework for nonlinear dimensionality reduction. Science, 290:2319–2323, Dec. 2000.
 17 A. Zomorodian and G. Carlsson. Computing persistent homology. Discrete and Computational Geometry, 33(2):249–274, 2005.

