\documentclass[amscd, amssymb, verbatim]{amsart}[12pt]

\usepackage{fullpage}

\usepackage{amssymb}

\usepackage{latexsym}

\usepackage{pb-diagram}

\usepackage[mathscr]{euscript}

\usepackage{graphicx}

\usepackage{subfigure}

\usepackage{placeins}

\usepackage{wrapfig}

\usepackage{hyperref}

\usepackage{appendix}

\usepackage{multind}

\usepackage{natbib}

\bibpunct{[}{]}{;}{a}{}{,~}

\theoremstyle{remark}

\newtheorem{exercise}{Exercise}

\theoremstyle{remark}

\newtheorem{exerciseSol}{Exercise}

\theoremstyle{remark}

\newtheorem*{solution}{Solution}

\DeclareMathOperator{\VR}{VR}

\DeclareMathOperator{\W}{W}

\DeclareMathOperator{\LW}{LW}

\newcommand{\N}{\mathbb{N}}

\newcommand{\R}{\mathbb{R}}

\newcommand{\Z}{\mathbb{Z}}

\newcommand{\cyclo}{\mathrm{C}_8\mathrm{H}_{16}}

\begin{document}

\title{Javaplex Tutorial}

\author{Henry Adams}

\email[Henry~Adams]{adams@math.colostate.edu}

\author{Andrew Tausz}

\email[Andrew~Tausz]{andrew.tausz@gmail.com}

\date{\today}

\maketitle

\tableofcontents

\setlength{\parindent}{0pt}

\setlength{\parskip}{8pt}

% TODO: Move dump example up to earlier section of tutorial?

% TODO: test tutorial after making changes to include examples finding representative cycles (For example, getDefaultSimplicialAlgorithm changed everywhere to getModularSimplicialAlgorithm).

% TODO: Add more figures with data points in blue and landmarks in red, such as the one in exercise_19, to this tutorial.

% TODO: the figures for the random trials in this tutorial don't currently match the numbers (number of simplices, etc).

%UPDATE: They do for witness_example.m and lazy_witness_example.m

% "Googlecode Issue": In the tutorial, explain the new feature that allows the user to specify the initial seed used for sequential maxmin landmark selection. This feature has been added as a result of Issue 23 (fixed).

% TODO: explain API class and how most commands for users are in this class (API does not contain under-the-hood implementation, but instead shows what the user can do w/ the library.) Check if Lazy Witness is in API; if not, perhaps it should be added.

% TODO: add example which is FLAT double torus. That is, select points at random from an octagon. Then, create explicit metric by identifying sides.

% UPDATE: added example for double torus embedded in R^3.

% TODO: is "ensureAllFaces" working as I expect? See representatives_example_HA and try changing dimensions.

% UPDATE: See Javaplex Issue 20 created by Andrew.

% TODO: add exercise numbers to tutorial wiki.

% TODO: Add zigzag examples to tutorial. See:

% http://code.google.com/p/javaplex/source/browse/#svn%2Ftrunk%2Fsrc%2Fmatlab%2Fexperimental%2Fvietoris_rips_bootstrap

% http://code.google.com/p/javaplex/source/browse/#svn%2Ftrunk%2Fsrc%2Fmatlab%2Fexperimental%2Fwitness_bootstrap

% http://arxiv.org/abs/1108.3545

% Also, see email from Anthony Bak 10-21-11

% TODO: ask Andrew: Move lazy witness to API so that the syntax is more similar to Vietoris-Rips and witness filtrations?

\section{Introduction}

\subsection{Javaplex}

Javaplex is a Java software package for computing the persistent homology of filtered simplicial complexes (or more generally, filtered chain complexes), with special emphasis on applications arising in topological data analysis \citep{Javaplex}. The main author is Andrew Tausz. Javaplex is a re-write of the JPlex package, which was written by Harlan Sexton and Mikael Vejdemo--Johansson. The main motivation for the development of Javaplex was the need for a flexible platform that supported new directions of research in topological data analysis and computational persistent homology. The website for Javaplex is \url{http://appliedtopology.github.io/javaplex/}, the documentation overview is at \url{https://github.com/appliedtopology/javaplex/wiki/Overview}, and the javadoc tree for the library is at \url{http://appliedtopology.github.io/javaplex/doc/}.

This tutorial is written for those using Javaplex with Matlab. However, one can run Javaplex without Matlab; see \url{https://github.com/appliedtopology/javaplex/wiki/Interoperability}.

If you are interested in Javaplex, then you may also be interested in the software package Dionysus by Dmitriy Morozov (\url{http://www.mrzv.org/software/dionysus}) or the software package Persus by Vidit Nanda (\url{http://www.sas.upenn.edu/~vnanda/perseus/index.html}).

% Below link currently broken:

%Some of the exercises in this tutorial are borrowed from Vin de Silva's {\em Plexercises}, available at \url{http://comptop.stanford.edu/u/programs/Plexercises2.pdf}.

Please email Henry at \texttt{adams@math.colostate.edu} or Andrew at \texttt{andrew.tausz@gmail.com} if you have questions about this tutorial.

\subsection{License}

Javaplex is an open source software package under the Open BSD License. The source code can be found at \url{https://github.com/appliedtopology/javaplex}.

\subsection{Installation for Matlab}

Open Matlab and check which version of Java is being used. In this tutorial, the symbol \texttt{>>} precedes commands to enter into your Matlab window.

\begin{quote} \begin{verbatim}

>> version -java

ans = Java 1.5.0_13 with Apple Inc. Java Hotspot(TM) Client VM mixed mode, sharing

\end{verbatim} \end{quote}

Javaplex requires version number \texttt{1.5} or higher.

To install Javaplex for Matlab, go to the latest release at \url{https://github.com/appliedtopology/javaplex/releases/latest/}. Download the zip file containing the Matlab examples, which should be called something like \texttt{matlab-examples-4.2.3.zip}. Extract the zip file. The resulting folder should be called \texttt{matlab\_examples}; make sure you know the location of this folder and that it is not inside the zip file.

In Matlab, change Matlab's ``Current Folder" to the directory \texttt{matlab\_examples} that you just extracted from the zip file. Run the \texttt{load\_javaplex.m} file.

\begin{quote} \begin{verbatim}

>> load_javaplex

\end{verbatim} \end{quote}

Installation is complete. Confirm that Javaplex is working properly with the following command.

\begin{quote} \begin{verbatim}

>> api.Plex4.createExplicitSimplexStream()

ans = edu.stanford.math.plex4.streams.impl.ExplicitSimplexStream@513fd4

\end{verbatim} \end{quote}

Your output should be the same except for the last several characters. Each time upon starting a new Matlab session, you will need to run \texttt{load\_javaplex.m}.

\subsection{Accompanying files}

The folder \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples}{\texttt{tutorial\_examples}} contains Matlab scripts, such as \\ %Space

\href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/explicit_simplex_example.m}{\texttt{explicit\_simplex\_example.m}} or \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/house_example.m}{\texttt{house\_example.m}}, which list all of the commands in this tutorial. This means that you don't need to type in each command individually. The folder \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples}{\texttt{tutorial\_examples}} also contains the Matlab data files, such as \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/pointsRange.mat}{\texttt{pointsRange.mat}} or \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/pointsTorusGrid.mat}{\texttt{pointsTorusGrid.mat}}, which are used in this tutorial. The folder \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_solutions}{\texttt{tutorial\_solutions}} contains the solution scripts, such as \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_solutions/exercise_1.m}{\texttt{exercise\_1.m}}, for all of the tutorial exercises. See Appendix \ref{A:solutions} for exercise solutions.

\section{Math review}

Below is a brief math review. For more details, see \citet{Armstrong}, \citet{ComputationalTopology}, \citet{TopologicalPersistence}, \citet{Hatcher}, and \citet{ComputingPersistent}.

\subsection{Simplicial complexes}

An abstract simplicial complex is given by the following data.

\begin{itemize}

\item A set $Z$ of vertices or 0-simplices.

\item For each $k\geq 1$, a set of $k$-simplices $\sigma = [z_0,z_1,\ldots,z_k]$, where $z_i\in Z$.

\item Each $k$-simplex has $k+1$ faces obtained by deleting one of the vertices. The following membership property must be satisfied: if $\sigma$ is in the simplicial complex, then all faces of $\sigma$ must be in the simplicial complex.

\end{itemize}

We think of 0-simplices as vertices, 1-simplices as edges, 2-simplices as triangular faces, and 3-simplices as tetrahedrons.

\subsection{Homology}

Betti numbers help describe the homology of a simplicial complex $X$. The value $Betti_k$, where $k\in \N$, is equal to the rank of the $k$-th homology group of $X$. Roughly speaking, $Betti_k$ gives the number of $k$-dimensional holes. In particular, $Betti_0$ is the number of connected components. For instance, a $k$-dimensional sphere has all Betti numbers equal to zero except for $Betti_0 = Betti_k = 1$.

\subsection{Filtered simplicial complexes}

A filtration on a simplicial complex $X$ is a collection of subcomplexes $\{X(t) \ |\ t\in \R\}$ of $X$ such that $X(t) \subset X(t')$ whenever $t\leq t'$. The filtration value of a simplex $\sigma \in X$ is the smallest $t$ such that $\sigma \in X(t)$. In Javaplex, filtered simplicial complexes (or more generally filtered chain complexes) are called streams.

% Introduce filtered chain complexes?

\subsection{Persistent homology}

Betti intervals help describe how the homology of $X(t)$ changes with $t$. A $k$-dimensional Betti interval, with endpoints [$t_{start}$, $t_{end}$), corresponds roughly to a $k$-dimensional hole that appears at filtration value $t_{start}$, remains open for $t_{start} \leq t < t_{end}$, and closes at value $t_{end}$. We are often interested in Betti intervals that persist for a long filtration range.

Persistent homology depends heavily on functoriality: for $t\leq t'$, the inclusion $i:X(t)\to X(t')$ of simplicial complexes induces a map $i_*:H_k(X(t))\to H_k(X(t'))$ between homology groups.

\section{Explicit simplex streams}\label{explicitStream}

In Javaplex, filtered simplicial complexes (or more generally filtered chain complexes) are called streams. The class ExplicitSimplexStream allows one to build a simplicial complex from scratch. In Section \ref{sfpc} we will learn about other automated methods of generating simplicial complexes; namely the Vietoris--Rips, witness, and lazy witness constructions.

\subsection{Computing homology}

You should change your current Matlab directory to \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples}{\texttt{tutorial\_examples}}, perhaps using the following command.

\begin{quote} \begin{verbatim}

>> cd tutorial_examples

\end{verbatim} \end{quote}

The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/explicit_simplex_example.m}{\texttt{explicit\_simplex\_example.m}}, which is in the folder \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples}{\texttt{tutorial\_examples}}. You may copy and paste commands from this script into the Matlab window, or you may run the entire script at once with the following command.

\begin{quote} \begin{verbatim}

>> explicit_simplex_example

\end{verbatim} \end{quote}

{\em Circle example.} Let's build a simplicial complex homeomorphic to a circle. We have three 0-simplices: [0], [1], [2], and three 1-simplices: [0,1], [0,2], [1,2].

\begin{figure}[htp]

\begin{center}

\includegraphics[width=1in]{s1.pdf}

\end{center}

\end{figure}

\FloatBarrier

To build a simplicial complex in Javaplex we simply build a stream in which all filtration values are zero. First we create an empty explicit simplex stream. Many command lines in this tutorial will end with a semicolon to supress unwanted output.

\begin{quote} \begin{verbatim}

>> stream = api.Plex4.createExplicitSimplexStream();

\end{verbatim} \end{quote}

Next we add simplicies using the methods \texttt{addVertex} and \texttt{addElement}. The first creates a vertex with a specified index, and the second creates a $k$-simplex (for $k >0$) with the specified array of vertices. Since we don't specify any filtration values, by default all added simplices will have filtration value zero.

\begin{quote} \begin{verbatim}

>> stream.addVertex(0);

>> stream.addVertex(1);

>> stream.addVertex(2);

>> stream.addElement([0, 1]);

>> stream.addElement([0, 2]);

>> stream.addElement([1, 2]);

>> stream.finalizeStream();

\end{verbatim} \end{quote}

After we are done building the complex, calling the method \texttt{finalizeStream} is necessary before working with this complex!

We print the total number of simplices in the complex.

\begin{quote} \begin{verbatim}

>> num_simplices = stream.getSize()

num_simplices = 6

\end{verbatim} \end{quote}

We create an object that will compute the homology of our complex. The first input parameter 3 indicates that homology will be computed in dimensions 0, 1, and 2 --- that is, in all dimensions strictly less than 3. The second input 2 means that we will compute homology with $\Z_2$ coefficients, and this input can be any prime number.

\begin{quote} \begin{verbatim}

>> persistence = api.Plex4.getModularSimplicialAlgorithm(3, 2);

\end{verbatim} \end{quote}

We compute and print the intervals.

\begin{quote} \begin{verbatim}

>> circle_intervals = persistence.computeIntervals(stream)

circle_intervals =

Dimension: 1

[0.0, infinity)

Dimension: 0

[0.0, infinity)

\end{verbatim} \end{quote}

This gives us the expected Betti numbers $Betti_0=1$ and $Betti_1=1$.

The persistence algorithm computing intervals can also find a representative cycle for each interval. However, there is no guarantee that the produced representative will be geometrically nice.

\begin{quote} \begin{verbatim}

>> circle_intervals = persistence.computeAnnotatedIntervals(stream)

circle_intervals =

Dimension: 1

[0.0, infinity): [1,2] + [0,2] + [0,1]

Dimension: 0

[0.0, infinity): [0]

\end{verbatim} \end{quote}

A representative cycle generating the single 0-dimensional homology class is [0], and a representative cycle generating the single 1-dimensional homology class is [1,2] + [0,2] + [0,1].

{\em 9-sphere example.} Let's build a 9-sphere, which is homeomorphic to the boundary of a 10-simplex. First we add a single 10-simplex to an empty explicit simplex stream. The result is not a simplicial complex because it does not contain the faces of the 10-simplex. We add all faces using the method \texttt{ensureAllFaces}. Then, we remove the 10-simplex using the method \texttt{removeElementIfPresent}. What remains is the boundary of a 10-simplex, that is, a 9-sphere.

\begin{quote} \begin{verbatim}

>> dimension = 9;

>> stream = api.Plex4.createExplicitSimplexStream();

>> stream.addElement(0:(dimension + 1));

>> stream.ensureAllFaces();

>> stream.removeElementIfPresent(0:(dimension + 1));

>> stream.finalizeStream();

\end{verbatim} \end{quote}

In the above, the \texttt{finalizeStream} method is used to ensure that the stream has been fully constructed and is ready for consumption by a persistence algorithm. It should be called every time after you build an explicit simplex stream.

We print the total number of simplices in the complex.

\begin{quote} \begin{verbatim}

>> num_simplices = stream.getSize()

num_simplices = 2046

\end{verbatim} \end{quote}

We get the persistence algorithm

\begin{quote} \begin{verbatim}

persistence = api.Plex4.getModularSimplicialAlgorithm(dimension + 1, 2);

\end{verbatim} \end{quote}

and compute and print the intervals.

\begin{quote} \begin{verbatim}

>> n_sphere_intervals = persistence.computeIntervals(stream)

n_sphere_intervals =

Dimension: 9

[0.0, infinity)

Dimension: 0

[0.0, infinity)

\end{verbatim} \end{quote}

This gives us the expected Betti numbers $Betti_0=1$ and $Betti_9=1$.

Try computing a representative cycle for each barcode.

\begin{quote} \begin{verbatim}

>> n_sphere_intervals = persistence.computeAnnotatedIntervals(stream)

\end{verbatim} \end{quote}

We don't display the output from this command in the tutorial, because the representative 9-cycle is very long and contains all eleven 9-simplices.

See Appendix \ref{A:solutions} for exercise solutions.

\begin{exercise}

Build a simplicial complex homeomorphic to the torus. Compute its Betti numbers. {\em Hint: You will need at least 7 vertices} \citep[page 107]{Hatcher}{\em . We recommend using a $3\times 3$ grid of 9 vertices.}

\end{exercise}

\begin{exercise}

Build a simplicial complex homeomorphic to the Klein bottle. Check that it has the same Betti numbers as the torus over $\Z_2$ coefficients but different Betti numbers over $\Z_3$ coefficients.

\end{exercise}

\begin{exercise}

Build a simplicial complex homeomorphic to the projective plane. Find its Betti numbers over $\Z_2$ and $\Z_3$ coefficients.

\end{exercise}

\subsection{Computing persistent homology}\label{computingPersistentHomology}

Let's build a stream with nontrivial filtration values. If your filtration values are not all integers, then please see the remark at the end of this section.

{\em House example.} The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/house_example.m}{\texttt{house\_example.m}}.

\begin{wrapfigure}{r}{1in}

\begin{center}

\includegraphics[width=1in]{houseFig.pdf}

\end{center}

\end{wrapfigure}

We build a house, with the vertices and edges on the square appearing at value 0, with the top vertex appearing at value 1, with the roof edges appearing at values 2 and 3, and with the roof 2-simplex appearing at value 7.

\begin{quote} \begin{verbatim}

>> stream = api.Plex4.createExplicitSimplexStream();

>> stream.addVertex(1, 0);

>> stream.addVertex(2, 0);

>> stream.addVertex(3, 0);

>> stream.addVertex(4, 0);

>> stream.addVertex(5, 1);

>> stream.addElement([1, 2], 0);

>> stream.addElement([2, 3], 0);

>> stream.addElement([3, 4], 0);

>> stream.addElement([4, 1], 0);

>> stream.addElement([3, 5], 2);

>> stream.addElement([4, 5], 3);

>> stream.addElement([3, 4, 5], 7);

>> stream.finalizeStream();

\end{verbatim} \end{quote}

We get the persistence algorithm with $\Z_2$ coefficients

\begin{quote} \begin{verbatim}

>> persistence = api.Plex4.getModularSimplicialAlgorithm(3, 2);

\end{verbatim} \end{quote}

and compute the intervals.

\begin{quote} \begin{verbatim}

>> intervals = persistence.computeIntervals(stream)

intervals =

Dimension: 1

[3.0, 7.0)

[0.0, infinity)

Dimension: 0

[1.0, 2.0)

[0.0, infinity)

\end{verbatim} \end{quote}

There are four intervals. The first is a $Betti_1$ interval, starting at filtration value 3 and ending at 7, with a representative cycle formed by the three edges of the roof. This 1-dimensional hole forms when edge $[4,5]$ appears at filtration value 3 and closes when 2-simplex $[3,4,5]$ appears at filtration value 7.

We can also store the intervals as Matlab matrices.

\begin{quote} \begin{verbatim}

>> intervals_dim0 = edu.stanford.math.plex4.homology.barcodes.BarcodeUtility...

.getEndpoints(intervals, 0, 0)

intervals_dim0 =

0 Inf

1 2

>> intervals_dim1 = edu.stanford.math.plex4.homology.barcodes.BarcodeUtility...

.getEndpoints(intervals, 1, 0)

intervals_dim1 =

0 Inf

3 7

\end{verbatim} \end{quote}

The second input of this command is the dimension of the intervals, and the third input is a Boolean flag: 0 to include infinite intervals, and 1 to exclude infinite intervals.

We compute a representative cycle for each barcode.

\begin{quote} \begin{verbatim}

>> intervals = persistence.computeAnnotatedIntervals(stream)

intervals =

Dimension: 1

[3.0, 7.0): [4,5] + [3,4] + -[3,5]

[0.0, infinity): [1,4] + [2,3] + [1,2] + [3,4]

Dimension: 0

[1.0, 2.0): -[3] + [5]

[0.0, infinity): [1]

\end{verbatim} \end{quote}

One $Betti_0$ interval and one $Betti_1$ interval are semi-infinite.

\begin{quote} \begin{verbatim}

>> infinite_barcodes = intervals.getInfiniteIntervals()

infinite_barcodes =

Dimension: 1

[0.0, infinity): [1,4] + [2,3] + [1,2] + [3,4]

Dimension: 0

[0.0, infinity): [1]

\end{verbatim} \end{quote}

We can print the Betti numbers at the largest filtration value (7 in this case) as an array

\begin{quote} \begin{verbatim}

>> betti_numbers_array = infinite_barcodes.getBettiSequence()

betti_numbers_array =

1

1

\end{verbatim} \end{quote}

or as a list with entries of the form $k: Betti_k$.

\begin{quote} \begin{verbatim}

>> betti_numbers_string = infinite_barcodes.getBettiNumbers()

betti_numbers_string = {0: 1, 1: 1}

\end{verbatim} \end{quote}

The Matlab function \texttt{plot\_barcodes.m} lets us display the intervals as Betti barcodes. The Matlab structure array \texttt{options} contains different options for the plot. We choose the filename \texttt{house} and we choose the maximum filtration value for the plot to be eight.

\begin{quote} \begin{verbatim}

>> options.filename = 'house';

>> options.max_filtration_value = 8;

>> plot_barcodes(intervals, options);

\end{verbatim} \end{quote}

The file \texttt{house.png} is saved to your current directory.

\begin{figure}[htp]

\begin{center}

\includegraphics[width=6in]{house.png}

\end{center}

\end{figure}

\FloatBarrier

The filtration values are on the horizontal axis. The $Betti_k$ number of the stream at filtration value $t$ is the number of intervals in the dimension $k$ plot that intersect a vertical line through $t$. Check that the displayed intervals agree with the filtration values we built into the house stream. At value 0, a connected component and a 1-dimensional hole form. At value 1, a second connected component appears, which joins to the first at value 2. A second 1-dimensional hole forms at value 3, and closes at value 7.

{\em Remark.} The methods \texttt{addElement} and \texttt{removeElementIfPresent} do not necessarily enforce the definition of a stream. They allow us to build inconsistent complexes in which some simplex $\sigma \in X(t)$ contains a subsimplex $\sigma' \notin X(t)$, meaning that $X(t)$ is not a simplicial complex. The method \texttt{validateVerbose} returns \texttt{1} if our stream is consistent and returns \texttt{0} with explanation if not.

\begin{quote} \begin{verbatim}

>> stream.validateVerbose()

ans = 1

>> stream.addElement([1, 4, 5], 0);

>> stream.validateVerbose()

Filtration index of face [4,5] exceeds that of element [1,4,5] (3 > 0)

Stream does not contain face [1,5] of element [1,4,5]

ans = 0

\end{verbatim} \end{quote}

{\em Remark.} If you want to use filtration values that are not integers, then you first need to specify an upper bound on the filtration values in your complex. This is demonstrated below, where the non-integer filtration value is 17.23 and the upper bound is 100.

\begin{quote} \begin{verbatim}

>> stream = api.Plex4.createExplicitSimplexStream(100);

>> stream.addVertex(1, 17.23);

>> stream.finalizeStream();

\end{verbatim} \end{quote}

%In Section \ref{sfpc} we will create simplex streams that are not also explicit simplex streams. To display or edit such streams, we will first need to use the method \texttt{makeExplicit}. See Exercise \ref{ripsExpl}.

% Henry: At the moment there is no ``make Explicit'' function.

\section{Point cloud data}

A point cloud is a finite metric space, that is, a finite set of points equipped with a notion of distance. One can create a Euclidean metric space by specifying the coordinates of points in Euclidean space, or one can create an explicit metric space by specifying all pairwise distances between points. In Section \ref{sfpc} we will learn how to build streams from point cloud data.

\subsection{Euclidean metric spaces}\label{euc} The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/pointcloud_example.m}{\texttt{pointcloud\_example.m}}.

{\em House example.} Let's give Euclidean coordinates to the points of our house.

\vspace{-3mm}

\begin{figure}[htb]

\centering

\includegraphics[width=1.55in]{houseCoord.pdf}

\vspace{-3mm}

\caption{The house point cloud}

\label{fig:housePointCloud}

\end{figure}

\FloatBarrier

You can enter these coordinates manually.

\begin{quote} \begin{verbatim}

>> point_cloud = [-1,0; 1,0; 1,2; -1,2; 0,3]

point_cloud =

-1 0

1 0

1 2

-1 2

0 3

\end{verbatim} \end{quote}

Or alternatively, these coordinates are stored as a Javaplex example.

\begin{quote} \begin{verbatim}

>> point_cloud = examples.PointCloudExamples.getHouseExample();

\end{verbatim} \end{quote}

We create a metric space using these coordinates. The input to the \texttt{EuclideanMetricSpace} method is a matrix whose $i$-th row lists the coordinates of the $i$-th point.

\begin{quote} \begin{verbatim}

>> m_space = metric.impl.EuclideanMetricSpace(point_cloud);

\end{verbatim} \end{quote}

We can return the coordinates of a specific point. Note the points are indexed starting at 0.

\begin{quote} \begin{verbatim}

>> m_space.getPoint(0)

ans =

-1

0

>> m_space.getPoint(2)

ans =

1

2

\end{verbatim} \end{quote}

A metric space can return the distance between any two points.

\begin{quote} \begin{verbatim}

>> m_space.distance(m_space.getPoint(0), m_space.getPoint(2))

ans = 2.8284

\end{verbatim} \end{quote}

{\em Figure~8 example.} We select 1,000 points randomly from a figure eight, that is, the union of unit circles centered at $(0,1)$ and $(0,-1)$.

\begin{quote} \begin{verbatim}

>> point_cloud = examples.PointCloudExamples.getRandomFigure8Points(1000);

\end{verbatim} \end{quote}

We plot the points.

\begin{quote} \begin{verbatim}

>> figure

>> scatter(point_cloud(:,1), point_cloud(:,2), '.')

>> axis equal

\end{verbatim} \end{quote}

{\em Torus example.} We select 2,000 points randomly from a torus in $\R^3$ with inner radius 1 and outer radius 2. The first input is the number of points, the second input is the inner radius, and the third input is the outer radius

\begin{quote} \begin{verbatim}

>> point_cloud = examples.PointCloudExamples.getRandomTorusPoints(2000, 1, 2);

\end{verbatim} \end{quote}

We plot the points.

\begin{quote} \begin{verbatim}

>> figure

>> scatter3(point_cloud(:,1), point_cloud(:,2), point_cloud(:,3), '.')

>> axis equal

>> view(60,40)

\end{verbatim} \end{quote}

{\em Sphere product example.} We select 1,000 points randomly from the unit torus $S^1 \times S^1$ in $\R^4$. The first input is the number of points, the second input is the dimension of each sphere, and the third input is the number of sphere factors.

\begin{quote} \begin{verbatim}

>> point_cloud = examples.PointCloudExamples.getRandomSphereProductPoints(1000, 1, 2);

\end{verbatim} \end{quote}

Plotting the third and fourth coordinates of each point shows a circle $S^1$.

\begin{quote} \begin{verbatim}

>> figure

>> scatter(point_cloud(:,3), point_cloud(:,4), '.')

>> axis equal

\end{verbatim} \end{quote}

\subsection{Explicit metric spaces}

We can also create a metric space from a distance matrix using the method \texttt{ExplicitMetricSpace}. For a point cloud in Euclidean space, this method is generally less convenient than the command \texttt{EuclideanMetricSpace}. However, method \texttt{ExplicitMetricSpace} can be used for a point cloud in an arbitrary (perhaps non-Euclidean) metric space.

The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/explicit_metric_space_example.m}{\texttt{explicit\_metric\_space\_example.m}}.

{\em House example.} The matrix \texttt{distances} summarizes the metric for our house points in Figure~\ref{fig:housePointCloud}: entry $(i,j)$ is the distance from point $i$ to point $j$.

\begin{quote} \begin{verbatim}

>> distances = [0,2,sqrt(8),2,sqrt(10);

2,0,2,sqrt(8),sqrt(10);

sqrt(8),2,0,2,sqrt(2);

2,sqrt(8),2,0,sqrt(2);

sqrt(10),sqrt(10),sqrt(2),sqrt(2),0]

distances =

0 2.0000 2.8284 2.0000 3.1623

2.0000 0 2.0000 2.8284 3.1623

2.8284 2.0000 0 2.0000 1.4142

2.0000 2.8482 2.0000 0 1.4142

3.1623 3.1623 1.4142 1.4142 0

\end{verbatim} \end{quote}

We create a metric space from this distance matrix.

\begin{quote} \begin{verbatim}

>> m_space = metric.impl.ExplicitMetricSpace(distances);

\end{verbatim} \end{quote}

We return the distance between points 0 and 2.

\begin{quote} \begin{verbatim}

>> m_space.distance(0, 2)

ans = 2.8284

\end{verbatim} \end{quote}

% Be careful: the constructor \texttt{DistanceData()} will accept matrices that fail to be symmetric, square, or nonnegative, creating ``metrics'' that do not satisfy the mathematical definiton. The triangle inequality is similarly easy to ignore. However, asking for the distance between a point and itself will always return zero, even if you input a distance matrix with non-zero diagonal entries.

\begin{exercise}\label{Ex:flatTorus}

One way to produce a torus is to take a square $[0, 1] \times [0, 1]$ and then identify opposite sides. This is called the flat torus. More explicitly, the flat torus is the quotient space

$$([0, 1] \times [0, 1]) / \sim,$$

where $(0, y) \sim (1, y)$ for all $y \in [0, 1]$ and where $(x, 0) \sim (x, 1)$ for all $x \in [0, 1]$. The Euclidean metric on $[0, 1] \times [0, 1]$ induces a metric on the flat torus. For example, in the induced metric on the flat torus, the distance between $(0, \frac{1}{2})$ and $(1, \frac{1}{2})$ is zero, since these two points are identified. The distance between $(\frac{1}{10}, \frac{1}{2})$ and $(\frac{9}{10}, \frac{1}{2})$ is $\frac{2}{10}$, by passing through the point $(0, \frac{1}{2}) \sim (1, \frac{1}{2})$.

Write a Matlab script or function that selects 1,000 random points from the square $[0, 1] \times [0, 1]$ and then computes the 1,000 $\times$ 1,000 distance matrix for these points under the induced metric on the flat torus. Create an explicit metric space from this distance matrix.

This exercise is continued by Exercise \ref{Ex:flatTorusLazy}.

% TODO: Gard Spreeman at MRC gave me a good idea for the Matlab function FlatTorusDistanceMatrix.m: just move one point to the center of the square, and then find the distance in the square to the other point.

\end{exercise}

\begin{exercise}\label{Ex:flatKlein}

One way to produce a Klein bottle is to take a square $[0, 1] \times [0, 1]$ and then identify opposite edges, with the left and right sides identified with a twist. This is called the flat Klein bottle. More explicitly, the flat Klein bottle is the quotient space

$$([0, 1] \times [0, 1]) / \sim,$$

where $(0, y) \sim (1, 1 - y)$ for all $y \in [0, 1]$ and where $(x, 0) \sim (x, 1)$ for all $x \in [0, 1]$. The Euclidean metric on $[0, 1] \times [0, 1]$ induces a metric on the flat Klein bottle. For example, in the induced metric on the flat Klein bottle, the distance between $(0, \frac{4}{10})$ and $(1, \frac{6}{10})$ is zero, since these two points are identified. The distance between $(\frac{1}{10}, \frac{4}{10})$ and $(\frac{9}{10}, \frac{6}{10})$ is $\frac{2}{10}$, by passing through the point $(0, \frac{4}{10}) \sim (1, \frac{6}{10})$.

Write a Matlab script or function that selects 1,000 random points from the square $[0, 1] \times [0, 1]$ and then computes the 1,000 $\times$ 1,000 distance matrix for these points under the induced metric on the flat Klein bottle. Create an explicit metric space from this distance matrix.

This exercise is continued by Exercise \ref{Ex:flatKleinLazy}.

\end{exercise}

\begin{exercise}\label{Ex:quotProjPlane}

One way to produce a projective plane is to take the unit sphere $S^2 \subset \R^3$ and then identify antipodal points. More explicitly, the projective plane is the quotient space

$$S^2 / (x \sim -x).$$

The Euclidean metric on $S^2$ induces a metric on the projective plane.

Write a Matlab script or function that selects 1,000 random points from the unit sphere $S^2 \subset \R^3$ and then computes the 1,000 $\times$ 1,000 distance matrix for these points under the induced metric on the projective plane. Create an explicit metric space from this distance matrix.

This exercise is continued by Exercise \ref{Ex:quotProjPlaneLazy}.

\end{exercise}

\section{Streams from point cloud data}\label{sfpc}

In Section \ref{explicitStream} we built streams explicitly, or by hand. In this section we construct streams from a point cloud $Z$. We build Vietoris--Rips, witness, and lazy witness streams. See \citet{WitnessComplexes} for additional information.

The Vietoris--Rips, witness, and lazy witness streams all take three of the same inputs: the maximum dimension $d_{max}$ of any included simplex, the maximum filtration value $t_{max}$, and the number of divisions $N$. These inputs allow the user to limit the size of the constructed stream, for computational efficiency. No simplices above dimension $d_{max}$ are included. The persistent homology of the resulting stream can be calculated only up to dimension $d_{max} - 1$ since homology in dimension $d_{max}-1$ depends on the boundary matrix from $d_{max}$-simplices to $(d_{max}-1)$-simplices. Also, instead of computing filtered simplcial complex $X(t)$ for all $t \geq 0$, we only compute $X(t)$ for

$$t \in \Biggl\{ 0,\ \frac{t_{max}}{N-1},\ \frac{2t_{max}}{N-1},\ \frac{3t_{max}}{N-1},\ \dots,\ \frac{(N-2)t_{max}}{N-1},\ t_{max} \Biggr\}.$$

The number of divisions $N$ is an optional input. If this input parameter is not specified, then the default value $N = 20$ is used.

When working with a new dataset, don't choose $d_{max}$ and $t_{max}$ too large initially. First get a feel for how fast the simplicial complexes are growing, and then raise $d_{max}$ and $t_{max}$ nearer to the computational limits. If you ever choose $d_{max}$ or $t_{max}$ too large and Matlab seems to be running forever, pressing the \texttt{control} and \texttt{c} buttons simultaneously may halt the computation. See also the remark in Section~\ref{ss:java-heap-size}.

%The first author is currently working with Jan Segert on interactive visualizations of the Vietoris--Rips and Witness filtrations for the Wolfram Demonstrations Project. We have preliminary drafts of the demonstrations. Please email Henry if you'd like to check them out; in particular, the witness filtrations can be hard to visualize.

\subsection{Vietoris--Rips streams}\label{Vietoris--Rips}

Let $d(\ \cdot\ ,\ \cdot \ )$ denote the distance between two points in metric space $Z$. A natural stream to build is the Vietoris--Rips stream. The complex $\VR(Z,t)$ is defined as follows:

\begin{itemize}

\item{the vertex set is $Z$.}

\item{for vertices $a$ and $b$, edge $[ab]$ is included in $\VR(Z,t)$ if $d(a,b) \leq t$.}

\item{a higher dimensional simplex is included in $\VR(Z,t)$ if all of its edges are.}

\end{itemize}

Note that $\VR(Z,t) \subset \VR(Z,t')$ whenever $t\leq t'$, so the Vietoris--Rips stream is a filtered simplicial complex. Since a Vietoris--Rips complex is the maximal simplicial complex that can be built on top of its 1-skeleton, it is an example of a {\em clique complex} or a {\em flag complex}.

The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/rips_example.m}{\texttt{rips\_example.m}}.

{\em House example.} Let's build a Vietoris--Rips stream from the house point cloud in Section \ref{euc}, where the metric space is $Z=\{(-1,0),(1,0),(1,2),(-1,2),(0,3)\}$. Note this stream is different than the explicit house stream we built in Section \ref{computingPersistentHomology}.

\begin{quote} \begin{verbatim}

>> max_dimension = 3;

>> max_filtration_value = 4;

>> num_divisions = 100;

>> point_cloud = examples.PointCloudExamples.getHouseExample();

>> stream = api.Plex4.createVietorisRipsStream(point_cloud, max_dimension, ...

max_filtration_value, num_divisions);

\end{verbatim} \end{quote}

The ellipses in the command above should be omitted; they are included only to indicate that this command continues onto the next line.

The order of the inputs is \texttt{createVietorisRipsStream(}$Z,\ d_{max},\ t_{max},\ N$\texttt{)}. For a Vietoris--Rips stream, the parameter $t_{max}$ is the maximum possible edge length. Since $t_{max} = 4$ is greater than the diameter ($\sqrt{10}$) of our point cloud, all edges will eventually form.

Since $d_{max} = 3$ we can compute up to second dimensional persistent homology.

\begin{quote} \begin{verbatim}

>> persistence = api.Plex4.getModularSimplicialAlgorithm(max_dimension, 2);

>> intervals = persistence.computeIntervals(stream);

\end{verbatim} \end{quote}

We display the Betti intervals. Parameter \texttt{options.max\_filtration\_value} is the largest filtration value to be displayed; typically \texttt{options.max\_filtration\_value} is chosen to be \texttt{max\_filtration\_value}. Parameter \texttt{options.max\_dimension} is the largest persistent homology dimension to be displayed; typically \texttt{options.max\_dimension} is chosen to be \texttt{max\_dimension - 1} because in a stream with simplices computed up to dimension $d_{max}$ we can only compute persistent homology up to dimension $d_{max} - 1$.

\begin{quote} \begin{verbatim}

>> options.filename = 'ripsHouse';

>> options.max_filtration_value = max_filtration_value;

>> options.max_dimension = max_dimension - 1;

>> plot_barcodes(intervals, options);

\end{verbatim} \end{quote}

The file \texttt{ripsHouse.png} is saved to your current directory.

\begin{figure}[htp]

\begin{center}

\includegraphics[width=6in]{ripsHouse.png}

\end{center}

\end{figure}

\FloatBarrier

Check that these plots are consistent with the Vietoris--Rips definition: edges $[3,5]$ and $[4,5]$ appear at filtration value $t = \sqrt{2}$; the square appears at $t = 2$; the square closes at $t = \sqrt{8}$.

%\begin{exercise}

%Change \texttt{ripsHouse} into an explicit stream

%

%\begin{quote} \begin{verbatim}

%>> ripsExpl = Plex.makeExplicit(ripsHouse);

%\end{verbatim} \end{quote}

%

%Check that you can display and edit stream \texttt{ripsExpl} using the methods of Section \ref{explicitStream}.

%\end{exercise}

% Is the above exercise no longer necessary, or no longer possible?

% Answer: At the moment it isn't possible. Do you think it is meaningful?

% Is there a way to print the simplices in a stream (and perhaps their filtration values?

{\em Torus example.} Try the following sequence of commands. We start with 400 points from a $20 \times 20$ grid on the unit torus $S^1 \times S^1$ in $\R^4$, and add a small amount of noise to each point.

% If I select points randomly from the torus instead from the regular grid, then the Vietoris--Rips stream does not seem to be able to recover the correct Betti numbers.

We build the Vietoris--Rips stream.

\begin{quote} \begin{verbatim}

>> max_dimension = 3;

>> max_filtration_value = 0.9;

>> num_divisions = 100;

\end{verbatim} \end{quote}

Make sure you are in the directory \texttt{tutorial\_examples} (you may need to enter the command \texttt{cd tutorial\_examples}), and then load the file \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/pointsTorusGrid.mat}{\texttt{pointsTorusGrid.mat}}. The matrix \texttt{pointsTorusGrid} appears in your Matlab workspace.

\begin{quote} \begin{verbatim}

>> load pointsTorusGrid.mat

>> point_cloud = pointsTorusGrid;

>> size(point_cloud)

ans = 400 4 % 400 points in dimension 4

>> stream = api.Plex4.createVietorisRipsStream(point_cloud, max_dimension, ...

max_filtration_value, num_divisions);

>> num_simplices = stream.getSize()

num_simplices = 82479

>> persistence = api.Plex4.getModularSimplicialAlgorithm(max_dimension, 2);

>> intervals = persistence.computeIntervals(stream);

>> options.filename = 'ripsTorus';

>> options.max_filtration_value = max_filtration_value;

>> options.max_dimension = max_dimension - 1;

>> options.side_by_side = true;

>> plot_barcodes(intervals, options);

\end{verbatim} \end{quote}

Setting the parameter \texttt{options.side\_by\_side} equal to \texttt{true} makes it such that the Betti barcodes of different dimensions are plotted side by side instead of above and below each other. The file \texttt{ripsTorus.png} is saved to your current directory.

\begin{figure}[htp]

\begin{center}

\includegraphics[width=6in]{ripsTorus.png}

\end{center}

\end{figure}

\FloatBarrier

The diameter of this torus (before adding noise) is $\sqrt{8}$, so choosing $t_{max} = 0.9$ likely will not show all homological activity. However, the torus will be reasonably connected by this time. Note the semi-infinite intervals match the correct numbers $Betti_0 = 1$, $Betti_1 = 2$, $Betti_2 = 1$ for a torus.

\begin{quote} \begin{verbatim}

>> infinite_barcodes = intervals.getInfiniteIntervals();

>> betti_numbers_array = infinite_barcodes.getBettiSequence()

betti_numbers_array =

1

2

1

\end{verbatim} \end{quote}

This example makes it clear that the computed ``semi-infinite'' intervals do not necessarily persist until $t = \infty$: in a Vietoris--Rips stream, once $t$ is greater than the diameter of the point cloud, the Betti numbers for $\VR(Z,t)$ will be $ Betti_0 = 1$, $Betti_1 = Betti_2 = ... = 0$. The computed semi-infinite intervals are merely those that persist until $t = t_{max}$.

{\em Remark.} We can build Vietoris--Rips streams not only on top of Euclidean point clouds, but also on top of explicit metric spaces. For example, if \texttt{m\_space} were an explicit metric space, then we could call a command such as

\begin{quote} \begin{verbatim}

>> stream = api.Plex4.createVietorisRipsStream(m_space, max_dimension, ...

max_filtration_value, num_divisions);

\end{verbatim} \end{quote}

\begin{exercise}

Slowly increase the values for $t_{max}$, $d_{max}$ and note how quickly the size of the Vietoris--Rips stream and the time of computation grow. Either increasing $t_{max}$ from 0.9 to 1 or increasing $d_{max}$ from 3 to 4 roughly doubles the size of the Vietoris--Rips stream.

\end{exercise}

\begin{exercise}

Find a planar dataset $Z \subset \R^2$ and a filtration value $t$ such that $\VR(Z,t)$ has nonzero $Betti_2$. Build a Vietoris--Rips stream to confirm your answer.

\end{exercise}

\begin{exercise}

Find a planar dataset $Z \subset \R^2$ and a filtration value $t$ such that $\VR(Z,t)$ has nonzero $Betti_6$. When building a Vietoris--Rips stream to confirm your answer, don't forget to choose $d_{max} = 7$.

\end{exercise}

\subsection{Landmark selection}\label{lands}

For larger datasets, if we include every data point as a vertex, as in the Vietoris--Rips construction, our streams will quickly contain too many simplices for efficient computation. The witness stream and the lazy witness stream address this problem. In building these streams, we select a subset $L \subset Z$, called landmark points, as the only vertices. All data points in $Z$ help serve as witnesses for the inclusion of higher dimensional simplices.

There are two common methods for selecting landmark points. The first is to choose the landmarks $L$ randomly from point cloud $Z$. The second is a greedy inductive selection process called sequential maxmin. In sequential maxmin, the first landmark is picked randomly from $Z$. Inductively, if $L_{i-1}$ is the set of the first $i-1$ landmarks, then let the $i$-th landmark be the point of $Z$ which maximizes the function $z \mapsto d(z, L_{i-1})$, where $d(z, L_{i-1})$ is the distance between the point $z$ and the set $L_{i-1}$.

Landmarks chosen using sequential maxmin tend to cover the dataset and to be spread apart from each other. A disadvantage is that outlier points tend to be selected. However, outlier points are less of an issue if one first takes dense core subsets as in Appendix \ref{A:core}. Sequential maxmin landmarks are used by \citet{Range} and \citet{KleinBottle}.

The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/landmark_example.m}{\texttt{landmark\_example.m}}.

{\em Figure~8 example.} We create a point cloud of 1,000 points from a figure eight.

\begin{quote} \begin{verbatim}

>> point_cloud = examples.PointCloudExamples.getRandomFigure8Points(1000);

\end{verbatim} \end{quote}

We create both a random landmark selector and a sequential maxmin landmark selector. These selectors will pick 100 landmarks each.

\begin{quote} \begin{verbatim}

>> num_landmark_points = 100;

>> random_selector = api.Plex4.createRandomSelector(point_cloud, num_landmark_points);

>> maxmin_selector = api.Plex4.createMaxMinSelector(point_cloud, num_landmark_points);

\end{verbatim} \end{quote}

We select 100 random landmarks and 100 landmarks via sequential maxmin. Note we need to increment the indices by 1 since Java uses 0-based arrays.

\begin{quote} \begin{verbatim}

>> random_points = point_cloud(random_selector.getLandmarkPoints() + 1, :);

>> maxmin_points = point_cloud(maxmin_selector.getLandmarkPoints() + 1, :);

\end{verbatim} \end{quote}

We plot the two sets of landmark points to see the difference between random and sequential maxmin landmark selection.

\begin{quote} \begin{verbatim}

>> subplot(1, 2, 1);

>> scatter(random_points(:,1), random_points(:, 2));

>> title('Random landmark selection');

>> subplot(1, 2, 2);

>> scatter(maxmin_points(:,1), maxmin_points(:, 2));

>> title('Maxmin landmark selection');

\end{verbatim} \end{quote}

\begin{figure}[htp]

\begin{center}

\includegraphics[width=5in]{fig8lands}

\end{center}

\end{figure}

\FloatBarrier

Sequential maxmin seems to do a better job of choosing landmarks that cover the figure eight and that are spread apart.

{\em Remark.} We can select landmark points not only from Euclidean point clouds but also from explicit metric spaces. For example, if \texttt{m\_space} is an explicit metric space, then we may select landmarks using a command such as the following.

\begin{quote} \begin{verbatim}

>> maxmin_selector = api.Plex4.createMaxMinSelector(m_space, num_landmark_points);

\end{verbatim} \end{quote}

Given point cloud $Z$ and landmark subset $L$, we define $\texttt{R} = \max_{z\in Z}\bigl\{d(z,L)\bigr\}$. Number \texttt{R} reflects how finely the landmarks cover the dataset. We often use it as a guide for selecting the maximum filtration value $t_{max}$ for a witness or lazy witness stream.

\begin{exercise}

Let $Z$ be the point cloud in Figure~\ref{fig:housePointCloud} from Section \ref{euc}, corresponding to the house point cloud. Suppose we are using sequential maxmin to select a set $L$ of 3 landmarks, and the first (randomly selected) landmark is $(1,0)$. Find by hand the other two landmarks in $L$.

\end{exercise}

\begin{exercise}

Let $Z$ be a point cloud and $L$ a landmark subset. Show that if $L$ is chosen via sequential maxmin, then for any $l_i,l_j\in L$, we have $d(l_i,l_j)\geq\texttt{R}$.

\end{exercise}

\subsection{Witness streams}

Suppose we are given a point cloud $Z$ and landmark subset $L$. Let $m_k(z)$ be the distance from a point $z \in Z$ to its $(k+1)$-th closest landmark point. The witness stream complex $\W(Z,L,t)$ is defined as follows.

\begin{itemize}

\item{the vertex set is $L$.}

\item{for $k>0$ and vertices $l_i$, the $k$-simplex $[l_0 l_1 ... l_k]$ is in $\W(Z,L,t)$ if all of its faces are, and if there exists a witness point $z \in Z$ such that $$\max\bigl\{d(l_0,z), d(l_1,z), ..., d(l_k,z)\bigr\} \leq t+m_k(z).$$ }

\end{itemize}

Note that $\W(Z,L,t) \subset \W(Z,L,t')$ whenever $t\leq t'$, so the witness stream is a filtered simplicial complex. Note that a landmark point can serve as a witness point.

\begin{exercise}\label{Ex:witnessHouse}

Let $Z$ be the point cloud in Figure~\ref{fig:housePointCloud} from Section \ref{euc}, corresponding to the house point cloud. Let $L = \{(1,0),(0,3),(-1,0)\}$ be the landmark subset. Find by hand the filtration value for the edge between vertices $(1,0)$ and $(0,3)$. Which point or points witness this edge? What is the filtration value for the lone 2-simplex $[(1,0),(0,3),(-1,0)]$?

\end{exercise}

The Matlab script corresponding to this section is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/witness_example.m}{\texttt{witness\_example.m}}.

{\em Torus example.} Let's build a witness stream instance for 10,000 random points from the unit torus $S^1 \times S^1$ in $\R^4$, with 50 sequential maxmin landmarks.

\begin{quote} \begin{verbatim}

>> num_points = 10000;

>> num_landmark_points = 50;

>> max_dimension = 3;

>> num_divisions = 100;

>> point_cloud = examples.PointCloudExamples.getRandomSphereProductPoints(num_points, ...

1, 2);

>> landmark_selector = api.Plex4.createMaxMinSelector(point_cloud, num_landmark_points);

\end{verbatim} \end{quote}

The next command returns the landmark covering measure \texttt{R} from Section \ref{lands}. Often the value for $t_{max}$ is chosen in proportion to \texttt{R}.

\begin{quote} \begin{verbatim}

>> R = landmark_selector.getMaxDistanceFromPointsToLandmarks()

R = 0.7033 % Generally close to 0.7

>> max_filtration_value = R / 8;

\end{verbatim} \end{quote}

We create the witness stream.

\begin{quote} \begin{verbatim}

>> stream = api.Plex4.createWitnessStream(landmark_selector, max_dimension, ...

max_filtration_value, num_divisions);

>> num_simplices = stream.getSize()

num_simplices = 1164 % Generally close to 1200

\end{verbatim} \end{quote}

We plot the Betti intervals.

\begin{quote} \begin{verbatim}

>> persistence = api.Plex4.getModularSimplicialAlgorithm(max_dimension, 2);

>> intervals = persistence.computeIntervals(stream);

>> options.filename = 'witnessTorus';

>> options.max_filtration_value = max_filtration_value;

>> options.max_dimension = max_dimension - 1;

>> plot_barcodes(intervals, options);

\end{verbatim} \end{quote}

The file \texttt{witnessTorus.png} is saved to your current directory.

\begin{figure}[htp]

\begin{center}

\includegraphics[width=6in]{witnessTorus.png}

\end{center}

\end{figure}

\FloatBarrier

The idea of persistent homology is that long intervals should correspond to real topological features, whereas short intervals are considered to be noise. The plot above shows that for a long range, the torus numbers $Betti_0 = 1$, $Betti_1 = 2$, $Betti_2 = 1$ are obtained. Your plot should contain a similar range.

The witness stream above contains approximately 2,000 simplices, fewer than the approximately 80,000 simplices in the Vietoris--Rips stream from the torus example in Section \ref{Vietoris--Rips}. This is despite the fact that we started with a point cloud of 100,000 points in the witness case, but of only 400 points in the Vietoris--Rips case. This supports our belief that the witness stream returns good results at lower computational expense.

\subsection{Lazy witness streams}

A lazy witness stream is similar to a witness stream. However, there is an extra parameter $\nu$, typically chosen to be 0, 1, or 2, which helps determine how the lazy witness complexes $\LW_\nu(Z,L,t)$ are constructed. See \citet{WitnessComplexes} for more information.

Suppose we are given a point cloud $Z$, landmark subset $L$, and parameter $\nu\in\N$. If $\nu = 0$, let $m(z) = 0$ for all $z\in Z$. If $\nu >0$, let $m(z)$ be the distance from $z$ to the $\nu$-th closest landmark point. The lazy witness complex $\LW_\nu(Z,L,t)$ is defined as follows.

\begin{itemize}

\item{the vertex set is $L$.}

\item{for vertices $a$ and $b$, edge $[ab]$ is in $\LW_\nu(Z,L,t)$ if there exists a witness $z \in Z$ such that $$\max\bigl\{d(a,z), d(b,z)\bigr\} \leq t+m(z).$$}

\item{a higher dimensional simplex is in $\LW_\nu(Z,L,t)$ if all of its edges are.}

\end{itemize}

Note that $\LW_\nu(Z,L,t) \subset \LW_\nu(Z,L,t')$ whenever $t\leq t'$, so the lazy witness stream is a filtered simplicial complex. The adjective {\em lazy} refers to the fact that the lazy witness complex is a flag complex: since the 1-skeleton determines all higher dimensional simplices, less computation is involved.

\begin{exercise}

Let $Z$ be the point cloud in Figure~\ref{fig:housePointCloud} from Section \ref{euc}, corresponding to the house point cloud. Let $L = \{(1,0),(0,3),(-1,0)\}$ be the landmark subset. Let $\nu = 1$. Find by hand the filtration value for the edge between vertices $(1,0)$ and $(0,3)$. Which point or points witness this edge? What is the filtration value for the lone 2-simplex $[(1,0),(0,3),(-1,0)]$?

\end{exercise}

\begin{exercise}

Repeat the above exercise with $\nu = 0$ and with $\nu = 2$.

\end{exercise}

\begin{exercise}

Check that the 1-skeleton of a witness complex $\W(Z,L,t)$ is the same as the 1-skeleton of a lazy witness complex $\LW_2(Z,L,t)$. As a consequence, $\LW_2(Z,L,t)$ is the flag complex of $\W(Z,L,t)$.

\end{exercise}

{\em 2-sphere example.} The Matlab script corresponding to this example is \href{https://github.com/appliedtopology/javaplex/tree/master/src/matlab/for_distribution/tutorial_examples/lazy_witness_example.m}{\texttt{lazy\_witness\_example.m}}.

We use parameter $\nu = 1$.

\begin{quote} \begin{verbatim}

>> max_dimension = 3;

>> num_points = 1000;

>> num_landmark_points = 50;

>> nu = 1;

>> num_divisions = 100;

>> point_cloud = examples.PointCloudExamples.getRandomSpherePoints(num_points, ...

max_dimension - 1);

>> landmark_selector = api.Plex4.createMaxMinSelector(point_cloud, num_landmark_points);

\end{verbatim} \end{quote}