text:
\begin{document}

	\def\spacingset#1{\renewcommand{\baselinestretch}
		{#1}\small\normalsize} \spacingset{1}

	\if0\blind
	{
		\title{\bf Hypothesis Testing for Shapes using Vectorized Persistence Diagrams}
		\author[1]{Chul Moon\thanks{Corresponding author. Email: chulm@smu.edu}}
		\author[2]{Nicole A. Lazar}

		\affil[1]{Southern Methodist University}
		\affil[2]{Pennsylvania State University}

		\date{}

		\maketitle
	} \fi

	\if1\blind
	{
		\bigskip
		\bigskip
		\bigskip
		\begin{center}
			{\LARGE\bf Title}
		\end{center}
		\medskip
	} \fi

	\bigskip
	\begin{abstract}
		Topological data analysis involves the statistical characterization of the shape of data. Persistent homology is a primary tool of topological data analysis, which can be used to analyze those topological
		features and perform statistical inference.
		In this paper, we present a two-stage hypothesis test for vectorized persistence diagrams.
		The first stage filters elements in the vectorized persistence diagrams to reduce false positives.
		The second stage consists of multiple hypothesis tests, with false positives controlled by false discovery rates.
		We demonstrate applications of the proposed procedure on simulated point clouds and three-dimensional rock image data.
		Our results show that the proposed hypothesis tests can provide flexible and informative inferences on the shape of data with lower computational cost compared to the permutation test.

	\end{abstract}

	\noindent
	{\it Keywords: Topological Data Analysis, Persistent Homology, Statistical Inference, Hypothesis Testing}
	\vfill

	\newpage
	\spacingset{1.45}

	\section{Introduction}

	Modern science is facing a rapid increase in the volume of data as well as their complexity. Non-standard data types such as functional \citep{Ramsey2005}, manifold \citep{Genovese2012}, and object-oriented \citep{Marron2014} data have become more common.
	New methodologies need to be developed to analyze and gain useful information from these data.
	For example, Figure~\ref{fig:rocks} shows Micro-CT images obtained from two different types of sandstones.
	The properties of rocks, such as permeability, play an important role in earth science fields, including hydrogeology and petroleum engineering.
	However, such properties are difficult to estimate; one of the reasons is that they are influenced by parameters linked to pore geometry \citep{Hommel2018}.

	Because it is not easy to quantify pore geometry, various indirect measures are used instead, such as porosity, tortuosity, and specific surface area \citep{Bernabe1998,XIONG2016}.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/bent.png}
			\caption{Bentheimer}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/dodd.png}
			\caption{Doddington}
		\end{subfigure}

		\caption{Micro-CT images of Bentheimer and Doddington sandstones  \citep{icl}.}
		\label{fig:rocks}
	\end{figure}

	Topological data analysis (TDA) is a recent addition to the analytical toolbox, that quantifies the shape of data using their topological features, such as connected components and holes.
	A primary TDA tool is persistent homology, which analyzes topological features of data in various scales \citep{Edelsbrunner2008,Carlsson2009}.
	Persistent homology provides a numeric descriptor of the shape of data that is robust to noise and insensitive to metrics \citep{Chazal2017}.
	It has been applied to a wide variety of data sets, including atomic configurations \citep{Hiraoka2016}, biomolecular structures \citep{Cang2017}, brain arteries \citep{Bendich2016}, and porous materials \citep{Robins2016}.

	Using results from TDA, it is possible to move toward statistical inference on the shape of data.

	For example, by computing persistent homology of the two rocks in Figure~\ref{fig:rocks}, we can obtain numeric outputs that describe their shapes, such as size, connectivity, and structure of pores.

	We can then, in principle, use this output to conduct hypothesis tests to distinguish rock types based on their shapes such as pore shapes and connectivities.

	Most statistical and machine learning methods, however, cannot be directly applied to persistent homology results. First, numeric outputs of persistent homology are not vectors. The numeric values that we obtain are sets of intervals that describe how topological features persist. Also, the number of intervals produced differs from dataset to dataset. Therefore, we cannot use approaches developed for symbolic data \citep{Billard2006}.

	Various approaches have been suggested to represent persistent homology results in different spaces while preserving the summarized topological information; in Euclidean space \citep{Bendich2016tracking, Adams2017},
	reproducing kernel Hilbert space (RKHS) \citep{Reininghaus2015,Kusano2016},
	and $L^2$-space \citep{Bubenik2015}.

	Several methods have been proposed to make statistical inference using persistent homology results including
	confidence interval using bootstrap \citep{Fasy2014},
	Bayesian approaches \citep{Maroulas2019},
	linear models \citep{Obayashi2018},
	and numerical summaries using polynomials \citep{Adcock2016} and tropical geometry \citep{Kalisnik2018}.
	There are a few authors who describe hypothesis testing procedures for persistent homology; permutation tests using a pairwise distance of persistence diagrams \citep{Robinson2017,cericola2018} and

	functional summaries \citep{Chen2015,Berry2020}, multiple hypothesis testing procedures using uniform point cloud data \citep{Vejdemo2018}, tests using mean persistence landscapes \citep{Bubenik2015} and kernel hypothesis tests in RKHS \citep{Kusano2019}.

	In this paper, we propose a two-stage hypothesis test of filtering and testing for persistent homology features represented as vectors in Euclidean space.

	In the filtering stage, elements that may generate false positives are dropped.
	In the testing stage, we conduct multiple hypothesis tests and control the false discovery rate (FDR).

	Two-stage hypothesis tests in our paper are implemented to reduce false positives whereas they have been used to enhance the detection power for high-dimensional data, including microarray \citep{Hackstadt2009,Tritchler2009} and genome-wide associations \citep{Murcray2008,Kooperberg2008}.

	The rest of the paper is organized as follows.
	In Section 2, we provide background on persistent homology and issues related to hypothesis tests for vectorized persistence diagrams.
	Section 3 introduces the proposed two-stage hypothesis test procedures.
	Section 4 presents the hypothesis testing results for simulated point cloud data and rock image data.
	Finally, in Section 5, we discuss our main contributions and future directions.

	\section{Topological Data Analysis Background}

	\subsection{Homology and Persistent Homology}

	Consider a shape of interest that is given as a topological space $\mathcal{X}$. Homology analyzes $\mathcal{X}$ by examining its $k$-dimensional holes: connected components (zero-dimensional holes), loops (one-dimensional holes), and voids (two-dimensional holes).
	These $k$-dimensional holes are often called the topological features of the space.
	The $k$-th homology group $H_k(\mathcal{X})$ describes the $k$-dimensional holes in $\mathcal{X}$. The rank of $H_k(\mathcal{X})$, the Betti number $\beta_k$, counts the number of $k$-dimensional holes.
	Figure~\ref{fig:homology} illustrates three spaces and their corresponding Betti numbers.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}
			\filldraw (0,0) circle (2.5pt);
			\filldraw (0.6,1) circle (1.5pt);
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}
			\draw (2.5,0) circle(0.6);
			\draw (1.5,0) -- (0,0) -- (0,1.5) -- (1.5,1.5) -- cycle;
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}[scale=0.7]
			\draw (-1,0) to[bend left] (1,0);
			\draw (-1.2,.1) to[bend right] (1.2,.1);
			\draw[rotate=0] (0,0) ellipse (80pt and 50pt);
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0,-1) {$\beta_0=2$};
			\node[text width=1.5cm] at (0,-1.5) {$\beta_1=0$};
			\node[text width=1.5cm] at (0,-2) {$\beta_2=0$};
			\end{tikzpicture}
			\caption{Two points}
		\end{subfigure}
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0,-1) {$\beta_0=2$};
			\node[text width=1.5cm] at (0,-1.5) {$\beta_1=2$};
			\node[text width=1.5cm] at (0,-2) {$\beta_2=0$};
			\end{tikzpicture}
			\caption{Square and circle}
		\end{subfigure}
		\begin{subfigure}{0.3\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0,-1) {$\beta_0=1$};
			\node[text width=1.5cm] at (0,-1.5) {$\beta_1=2$};
			\node[text width=1.5cm] at (0,-2) {$\beta_2=1$};
			\end{tikzpicture}
			\caption{Torus}
		\end{subfigure}
		\caption{Betti numbers of three topological spaces. $\beta_0$, $\beta_1$, and $\beta_2$ are the number of zero, one, and two dimensional holes in a space. For example, a torus has one connected component, two loops, and one void.}
		\label{fig:homology}
	\end{figure}

	Persistent homology keeps track of the homology of a sequence of topological spaces. Let $\emptyset \subseteq \mathcal{X}_0 \subseteq \mathcal{X}_1 \subseteq \cdots \subseteq \mathcal{X}_n =\mathcal{X}$ be a non-decreasing sequence of topological spaces. For example, Figure~\ref{fig:rips} shows a sequence of topological spaces called the Rips complexes $\mathcal{R}_{\epsilon}$. A Rips complex $\mathcal{R}_{\epsilon}$ is a set consisting of $k$-dimensional simplices, such as vertex, edge, and triangle, where the $k$-dimensional simplices correspond to $k+1$ points whose pairwise distances are smaller than $\epsilon$ \citep{Hatcher2002}.
	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.16\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0.5,-0.7) {$\mathcal{R}_{0.5}$};
			\filldraw (0,0) circle (1.5pt);
			\filldraw (1,0) circle (1.5pt);
			\filldraw (1.3,1) circle (1.5pt);
			\filldraw (0.4,1.3) circle (1.5pt);
			\filldraw (-1.05,1) circle (1.5pt);
			\end{tikzpicture}
		\end{subfigure}
		$\subseteq$
		\begin{subfigure}{0.16\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0.5,-0.7) {$\mathcal{R}_1$};
			\filldraw (0,0) circle (1.5pt);
			\filldraw (1,0) circle (1.5pt);
			\filldraw (1.3,1) circle (1.5pt);
			\filldraw (0.4,1.3) circle (1.5pt);
			\filldraw (-1.05,1) circle (1.5pt);
			\draw (0,0) -- (1,0);
			\draw (1.3,1) -- (0.4,1.3);
			\end{tikzpicture}
		\end{subfigure}
		$\subseteq$
		\begin{subfigure}{0.16\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0.5,-0.7) {$\mathcal{R}_{1.2}$};

			\filldraw (0,0) circle (1.5pt);
			\filldraw (1,0) circle (1.5pt);
			\filldraw (1.3,1) circle (1.5pt);
			\filldraw (0.4,1.3) circle (1.5pt);
			\filldraw (-1.05,1) circle (1.5pt);
			\draw (0,0) -- (1,0);
			\draw (1.3,1) -- (0.4,1.3);
			\draw (1,0) -- (1.3,1);
			\end{tikzpicture}
		\end{subfigure}
		$\subseteq$
		\begin{subfigure}{0.16\linewidth}
			\centering
			\begin{tikzpicture}
			\node[text width=1.5cm] at (0.5,-0.7) {$\mathcal{R}_{1.4}$};

			\filldraw (0,0) circle (1.5pt);
			\filldraw (1,0) circle (1.5pt);
			\filldraw (1.3,1) circle (1.5pt);
			\filldraw (0.4,1.3) circle (1.5pt);
			\filldraw (-1.05,1) circle (1.5pt);
			\draw (0,0) -- (1,0);
			\draw (1.3,1) -- (0.4,1.3);
			\draw (1,0) -- (1.3,1);
			\draw (0,0) -- (0.4,1.3);

			\end{tikzpicture}
		\end{subfigure}
			$\subseteq$
		\begin{subfigure}{0.16\linewidth}
			\centering
			\begin{tikzpicture}
		\node[text width=1.5cm] at (0.5,-0.7) {$\mathcal{R}_{1.5}$};
		\draw[fill=yellow]  (0,0) -- (1,0) -- (0.4,1.3) -- (0,0) -- cycle;
		\draw[fill=green]  (1,0) -- (0.4,1.3) -- (1.3,1) -- cycle;
		\filldraw (0,0) circle (1.5pt);
		\filldraw (1,0) circle (1.5pt);
		\filldraw (1.3,1) circle (1.5pt);
		\filldraw (0.4,1.3) circle (1.5pt);
		\filldraw (-1.05,1) circle (1.5pt);
		\draw (0,0) -- (1,0);
		\draw (1.3,1) -- (0.4,1.3);
		\draw (1,0) -- (1.3,1);
		\draw (0,0) -- (0.4,1.3);
		\draw (0,0) -- (-1.05,1);
			\end{tikzpicture}
		\end{subfigure}
		\caption{Sequence of Rips complexes. The dots, lines, and colored regions are vertices, edges, and triangles, respectively.}
		\label{fig:rips}
	\end{figure}

	The inclusion of $\mathcal{X}_i$ in $\mathcal{X}_j$ for $0 \leq i \leq j \leq n$ induces a map called the homomorphism between the homology groups, $\iota_k(i,j):H_k(\mathcal{X}_i) \rightarrow H_k(\mathcal{X}_j)$. Therefore, for a given sequence of topological spaces, we have a sequence of homomorphisms
	\begin{equation*}
	0=H(\mathcal{X}_0) \rightarrow H(\mathcal{X}_1) \rightarrow \cdots \rightarrow H(\mathcal{X}_n) = H(\mathcal{X}).
	\end{equation*}

	The images of homomorphisms include information about the birth and death of topological features in the sequence of topological spaces.

	The image of $\iota_k(i,j)$ includes information of $k$-dimensional topological features that appear before $i$ and disappear after $j$. For example, the loop in Figure~\ref{fig:rips} appears at $\mathcal{R}_{1.4}$ ($\epsilon=1.4$) and disappears at $\mathcal{R}_{1.5}$ ($\epsilon=1.5$). Information of the loop is included in the image of the dimension one homology group $\iota_1(i,j)$ where $1.4<i<j<1.5$.

	We define the $k$-dimensional $(i,j)$ persistent homology groups as the image of the homomorphism $\iota_k(i,j):H_k(\mathcal{X}_i) \rightarrow H_k(\mathcal{X}_j)$.

	The rank of the $k$-dimensional $(i,j)$ persistent homology groups is equal to the number of $k$-dimensional topological features at $\mathcal{X}_j$ that are born at or before $\mathcal{X}_i$.
	Here, the parameter that controls the sequential changes, such as $i$ or $j$, is called the filtration.

	From persistent homology, we can observe dynamics of when a specific $k$-dimensional hole appears (birth) and disappears (death) over the filtration.
	As a result, persistent homology is given as multiple intervals of (birth, death).

	A persistence diagram is a graphical representation of persistent homology. The (birth, death) intervals are plotted as points in $\mathbb{R}^2$ of birth ($x$-axis) and death ($y$-axis). Because the death of a topological feature comes after its birth, points are plotted above the 45-degree line. Figure~\ref{fig:pd} shows the persistence diagram of the Rips complexes in Figure~\ref{fig:rips} for dimensions 0 (connected components) and 1 (loops). For example, the blue dots in Figure~\ref{fig:pd} summarize $\epsilon$ values when the data points are connected by edges and the red triangle at $(1.4,1.5)$ in Figure~\ref{fig:pd} represents the loop created at $\epsilon=1.4$, which dies at $\epsilon=1.5$.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.5\textwidth}
			\begin{tikzpicture}[scale=0.9]
			\begin{axis}[
			axis lines = left,
			ymin=-0.05,
			ymax=2.05,
			xmin=-0.05,
			xmax=2.05,
			xlabel = Birth,
			ylabel = Death,
			xtick={0,0.4,...,2},
			ytick={0,0.4,...,2},
			legend pos=south east,
			scatter/classes={
				z={mark=*,draw=blue},
			    o={mark=triangle*,red}}
		    ]
			\addplot[scatter,only marks,
			scatter src=explicit symbolic,
			visualization depends on={\thisrow{myvalue} \as \myvalue}]
			table[meta=label] {
				x     y      label myvalue
				0  0.948 z 3
				0  1   z 3
				0  1.044 z 3
				0 1.486 z 3
				0 2 z 3
				1.36 1.43 o 3
			};
			\draw (0, 0) -- (190, 190);
		    \legend{Dim 0, Dim 1}
			\end{axis}
			\end{tikzpicture}
		\end{subfigure}
		\caption{Dimension 0 and 1 persistence diagram of Rips complexes in Figure~\ref{fig:rips}}\label{fig:pd}
	\end{figure}

	Persistent homology can explore different aspects of shapes when it is applied to different types of data. For example, when persistent homology is computed for the binary rock images, it can describe geometric properties or materials such as sizes of pores and grains, their connectivity, and pore-throat radius \citep{Robins2016}.
	Figure~\ref{fig:interpret0} and \ref{fig:interpret1} in Supplementary Material show examples of topological features of two-dimensional binary rock images and where such features appear in dimension 0 and 1 persistence diagrams.

	For more detailed explanation of persistent homology, see \cite{Hatcher2002}, \cite{Edelsbrunner2008}, and \cite{Ghrist2008}.

	\subsection{Vectorized Persistence Diagrams in Euclidean Space}

	Although persistence diagrams include information about topological features, it is difficult to use them as input in data analysis.
	This is because persistence diagrams are not numeric vectors, which most classical statistical and machine learning methods require.

	Ways to represent persistence diagrams as vectors include binning \citep{Bendich2016tracking} and persistence images \citep{Adams2017} in Euclidean space, persistence scale-space kernel \citep{Reininghaus2015} and persistence weighted Gaussian kernel \citep{Kusano2016} in RKHS, rank function \citep{Robins2016rank} and persistence landscapes \citep{Bubenik2015} in $L^2$ space.

	The vector representation of persistence diagrams in Euclidean space has some advantages. First, the representation in Euclidean space makes it easier to conduct prediction and classification by applying classical statistical models and machine learning methods such as regression and support vector machine \citep{Adams2017, Obayashi2018}. Second, it has a natural graphical summary unlike the other representations.

	\cite{Bendich2016tracking} suggest a binning method for a Euclidean vector representation. First, a persistence diagram $PD=\{(\text{birth},\text{death})\}$ is transformed into $P=\{(\text{birth},\text{death}-\text{birth})\}$. The transformed persistence diagrams are divided into a pre-determined number of bins. Then the number of topological features (points in persistence diagrams) included in each bin is counted and recorded as intensity values.
	However, the binning method is a stable vector representation of persistent diagrams; small changes in persistence diagrams may result in different vector representations.

	\cite{Adams2017} propose a different robust vectorization method called the persistence image. First, a persistence diagram $PD=\{(\text{birth},\text{death})\}$ is transformed into $P=\{(\text{birth},\text{death}-\text{birth})\}$, as before. With smoothing function $f_{(u,v)}$ for $(u,v) \in P$ and weight function $w$, the persistence surface $\rho_P$ of the transformed persistence diagram $P$ is
	\begin{equation*}
	\rho_{P} (x,y) = \sum\limits_{(u,v) \in P} f_{(u,v)}(x,y) \cdot w(u,v),
	\end{equation*}
	where $x$ and $y$ are the $(x,y)$-coordinates of the persistence surface.

	The persistence images are obtained by taking an integral of $\rho_{P}$ over a given grid. \cite{Adams2017} show that the persistence image is stable with respect to bottleneck and Wasserstein distances  \citep{Edelsbrunner2010}. In this paper, we use persistence images to represent persistence diagrams as vectors in Euclidean space.

	The smoothing and weight functions are selected to preserve information in the persistence diagrams. We use the Gaussian smoothing function $f_{(u,v)}(x,y | h) = \frac{1}{2\pi h^2 }\exp( -( (x-u)^2 + (y-v)^2 ) /2h^2)$ with $h=1.5*\text{(grid size)}$. One can choose different values of $h$, but the persistence image representation is robust to the choice of the variance of the Gaussian function \citep{Adams2017}.
	Table~\ref{tb:weight} summarizes the four weights we use in this paper. The linear weight assigns higher weights to persistent topological features whereas the constant weight gives equal weights to all features. The soft(hard) arctangent weight is closer to the constant(linear) weight. The arctangent weights are suggested in \cite{Kusano2016} for RKHS and extended to persistence images in \cite{Obayashi2018}.

	\begin{table}[!ht]
		\centering
		\begin{tabular}{|c|c|}
			\hline
			Name          & Weight                           \\ \hline
			Constant        & $w(u,v)=1$                         \\ \hline
			Soft arctangent & $w(u,v)=\arctan( 0.5(v)^{0.5})$ \\ \hline
			Hard arctangent & $w(u,v)=\arctan(v)$             \\ \hline
			Linear          & $w(u,v)=v$                         \\ \hline
		\end{tabular}
		\caption{Summary of the four weights. Here, $(u,v)$ are birth and death values in a transformed persistence diagram $P$.}\label{tb:weight}
	\end{table}

	Figure~\ref{fig:PI} illustrates the conversion steps from the persistence diagram to the persistence image. First, the (birth, death) pair is transformed to (birth, death$-$birth). Then the transformed persistence diagram is represented as the persistence image by using a smoothing function and a weight.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.25\textwidth}
			\begin{tikzpicture}[scale=0.6]
			\begin{axis}[
			axis lines = left,
			domain=0:1.05,
			xlabel = Birth,
			ylabel = Death,
			xtick={0,0.2,...,1},
			ytick={0,0.2,...,1},
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=0:1.05,
			samples=100,
			color=black,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				0.06  0.1    c
				0.02  0.17   c
				0.1   0.15   c
				0.12  0.15   c
				0.16  0.75   c

				0.42  0.73   c
				0.22  0.45   c
				0.45  0.57   c

				0.5   0.53    c
				0.53  0.75   c
				0.55  0.62   c

				0.73  0.85   c
				0.76  0.81    c
				0.85  0.92   c
				0.9   0.95    c
			};
			\end{axis}
			\end{tikzpicture}
		\end{subfigure}
		$\mbox{ }\Longrightarrow\mbox{ }$
		\begin{subfigure}{0.25\textwidth}

			\begin{tikzpicture}[scale=0.6]
			\begin{axis}[
			axis lines = left,
			domain=0:1.05,
			xlabel = Birth,
			ylabel = Death$-$Birth,
			xtick={0,0.2,...,1},
			ytick={0,0.2,...,1},
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=0:1.05,
			samples=100,
			color=white,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				0.06  0.04    c
				0.02  0.15   c
				0.1   0.05   c
				0.12  0.03   c
				0.16  0.59   c

				0.42  0.31   c
				0.22  0.23   c
				0.45  0.12   c

				0.5   0.03    c
				0.53  0.22   c
				0.55  0.07   c

				0.73  0.12   c
				0.76  0.05    c
				0.85  0.07   c
				0.9   0.05    c
			};
			\end{axis}
			\end{tikzpicture}
		\end{subfigure}
		$\mbox{ }\Longrightarrow\mbox{ }$
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/example}
		\end{subfigure}
		\caption{Steps of converting the persistence diagram into the persistence image. From the original persistence diagram (left) to the transformed persistence diagram (middle) to the persistence image (right). In this example, Gaussian smoothing function and soft arctangent weight are used for the persistence images.}\label{fig:PI}
	\end{figure}

	\subsection{Hypothesis Tests for Persistent Homology}

	One can consider a hypothesis test that compares the topological features of data using persistence diagrams.
	However, it is difficult to define probability distributions on the space of persistence diagrams because the space is infinite in dimension and has complicated geometry \citep{Robinson2017}.
	As a result, the hypothesis test methods for persistent homology have been studied in two ways: using permutations and represented persistence diagrams.

	First, the permutation-based tests are defined on the space induced by pairwise metrics of persistence diagrams and functional spaces.
	\cite{Robinson2017} propose a permutation test for a two-sample setting using persistence diagrams.
	Assume that we have persistence diagrams $PD_i$, $i\in\{1,2,\cdots,n\}$ obtained from two groups of size $n_1$ and $n_2$ ($n_1+n_2=n$).
	Let $G=\{I,J\}$ be labels for the two groups where $I,J \subset \{1, 2, \cdots, n \}$, $I \cap J = \emptyset$, $I \cup J = \{1, 2, \cdots, n \}$, $n(I)=n_1$, and $n(J)=n_2$.
	The permutation test uses a joint loss function $L$ of a label $G$ defined by
	\begin{equation*}
	L(G=\{I,J\})= \frac{2}{n_1(n_1-1)} \sum\limits_{i<j \in I} d(PD_i,PD_j) + \frac{2}{n_2(n_2-1)}\sum\limits_{i<j \in J} d(PD_i,PD_j),
	\end{equation*}
	where $d$ is a pairwise distance function for persistence diagrams such as bottleneck or Wasserstein distance. The joint loss value measures the similarity of persistence diagrams within groups. If persistence diagrams within groups are similar(different) to each other, then the joint loss value will be small(large).
	Algorithm~\ref{alg:rt} in Supplementary Material presents the permutation test procedure of \cite{Robinson2017}.
	\cite{cericola2018} extend the two-sample test scheme to multiple label group testing using the one-way analysis of variance (ANOVA) procedure and \cite{Vejdemo2018} propose procedures to control the multiple testing problem.
	Also, \cite{Chen2015} and \cite{Berry2020} propose permutation tests using functional representations of persistence diagrams.

	Second, hypothesis tests are conducted to the various persistence diagram representations.

	\cite{Bubenik2015} suggests various tests for persistence landscapes such as two-sample tests for mean persistence landscapes and Hotelling's $T^2$ test.
	\cite{Kusano2019} applies the two-sample kernel test of \cite{Gretton2006,Gretton2012} to persistence diagram representations using the persistent weighted Gaussian kernel.

	Existing hypothesis test approaches have some limitations. First, the computational costs for the permutation test using persistence diagrams are high. Computing pairwise distances for persistence diagrams is computationally expensive; $O(n_f^{1.5}\log n_f)$ for bottleneck distance \citep{Kerber2017} and $O(n_f^3)$ using the Hungarian algorithm \citep{Munkres1957} and $O(n_f^\alpha)$ using the auction algorithm \citep{Kerber2017} for Wasserstein distance, where $n_f$ is the number of features in the persistence diagrams and $\alpha$ is determined by approximation.
	Second, most hypothesis test methods except for the persistence landscape do not provide information on how each topological feature contributes to any observed differences.
	Therefore, it is difficult to identify which features play an essential role in the hypothesis test.

	Lastly, some methods do not provide an option for assigning weights to topological features that enables flexible interpretation. For example, the permutation test using persistence diagrams and the hypothesis test using persistence landscapes are only available under the fixed linear weight.

	Conducting a hypothesis test on persistence images can overcome the limitations; it can be performed with low computational cost, identify topological features that account for the differences, and implement weights.
	However, a naive application of hypothesis tests to persistence images can lead to an incorrect conclusion.

	\subsection{Naive Hypothesis Test for Persistence Images}
	\label{subsec:tstt}

	For the vectorized persistence diagrams in Euclidean space, we can measure the differences by comparing the mean values of vectors.
	We consider testing the difference between two groups of spaces ($\mathcal{X}_1$ and $\mathcal{X}_2$) for simplicity. It is straightforward to generalize the approach to more than two groups, similar to the transition from t-tests to ANOVA.
	Assume that $n_1$ and $n_2$ persistence diagrams are obtained from $\mathcal{X}_1$ and $\mathcal{X}_2$, respectively.
	Let $x_{(j,k)}^{i}$ be the $i$th element of group $j$'s $k$th persistence image, where  $i\in \{ 1, \cdots, m \}$, $j=1,2$, and $k\in \{ 1, \cdots, n_j \}$. Also, let $\mu_j^i=E_{k}(x^i_{(j,k)})$ be the mean of the $i$th element of group $j$'s persistence images.
	Then, the hypotheses for testing differences between two persistence images are
	\begin{eqnarray*}
		H_0:&& \mu_1^{i} = \mu_2^{i} \text{ for all } i \in \{ 1, \cdots, m\} \nonumber \\
		H_1:&& \mu_1^{i} \neq \mu_2^{i} \text{ for at least one } i \in \{ 1, \cdots, m\}.
		\label{eqn:hypo}
	\end{eqnarray*}

	However, it is not appropriate to implement these naive hypotheses directly. First, not all areas of the persistence images are equally interesting. There will most likely be multiple elements with a mean close to zero in both groups just because they are far from the diagonal or are in otherwise sparse regions.
	Second, the elements are not independent. The persistence image uses a smoothing function, so nearby elements can be dependent.

	To illustrate these issues, we reproduce a simulation study considered in \cite{Robinson2017}. We randomly sample 100 sets of 50 points from each of two different shapes: one circle of radius 1 (shape 1) and two circles of radius 0.9 and 1.1 (shape 2), as shown in Figure~\ref{fig:sim}.
	We compute persistent homology for 200 sets of point cloud data using the Rips complex.
	The difference between the two spaces is reflected in the one-dimensional persistence diagrams, given in the second from the right panel of Figure~\ref{fig:sim}. The dimension one interval for shape 1 ($ \text{Death}- \text{Birth} \approx1.3$) is longer than the longest dimension one interval for shape 2 ($ \text{Death}- \text{Birth} \approx1.1$). Also, shape 2 has several very short-lived dimension one intervals.
	The 200 one-dimensional persistence diagrams are converted to 200 persistence images.

		\begin{figure}[!ht]
		\centering

		\begin{subfigure}{0.2\textwidth}
			\centering
			\begin{tikzpicture}[y=.5cm, x=.5cm,font=\sffamily]
			\draw (0,2) circle (1.5cm);
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.24\textwidth}
			\centering
			\begin{tikzpicture}[scale=0.5]

			\begin{axis}[
			axis lines = left,
			domain=-1.2:1.2,
			xlabel = x,
			ylabel = y,
			axis equal,
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=-1.2:1.2,
			samples=100,
			color=white,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				0.7428219  0.66948913 c
				0.3385895  0.94093421 c
				-0.2571978 -0.96635877 c
				0.9766586  0.21479765 c
				-0.9878061  0.15568911 c
				-0.2739267  0.96175057 c
				-0.7914422 -0.61124407 c
				0.2394508 -0.97090850 c
				-0.8349805  0.55027957 c
				-0.2368550 -0.97154502 c
				0.5394487  0.84201850 c
				-0.5494803  0.83550667 c
				0.8730733  0.48758894 c
				0.5881371  0.80876126 c
				-0.8459769  0.53321957 c
				-0.4838730  0.87513821 c
				-0.7541640 -0.65668616 c
				-0.9157753  0.40169103 c
				-0.4527487 -0.89163819 c
				0.3751519  0.92696334 c
				0.9921755 -0.12485112 c
				-0.9994425  0.03338565 c
				0.4180952 -0.90840321 c
				-0.9997605 -0.02188412 c
				-0.5903125  0.80717481 c
				0.6363443  0.77140520 c
				0.7812096 -0.62426879 c
				-0.7736046 -0.63366856 c
				-0.8815971 -0.47200262 c
				0.6375824 -0.77038216 c
				0.4916157  0.87081229 c
				-0.9949560 -0.10031186 c
				0.9769235 -0.21358941 c
				-0.8313772 -0.55570852 c
				0.9203847  0.39101412 c
				0.5083239 -0.86116597 c
				-0.9142303 -0.40519487 c
				-0.1910957 -0.98157140 c
				0.8284818  0.56001596 c
				-0.9048752 -0.42567692 c
				0.1249782 -0.99215949 c
				-0.9731828 -0.23003293 c
				-0.9915786 -0.12950616 c
				0.2961873 -0.95512986 c
				-0.3613478  0.93243110 c
				0.9950132 -0.09974289 c
				0.2716518 -0.96239559 c
				-0.9805989 -0.19602506 c
				-0.2699407  0.96287695 c
				0.3747664  0.92711926 c
			};
			\end{axis}
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.24\textwidth}
			\centering
			\begin{tikzpicture}[scale=0.5]
			\begin{axis}[
			axis lines = left,
			domain=0:2.1,
			xlabel = Birth,
			ylabel = Death$-$Birth,
			xtick={0,0.4,...,2},
			ytick={0,0.4,...,2},
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=0:2.1,
			samples=100,
			color=white,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				0.3830768 1.354    c
			};
			\end{axis}
			\end{tikzpicture}
		\end{subfigure}
		\begin{subfigure}{0.29\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/exone_new}

			\label{subfig:ttesta}
		\end{subfigure}

		\begin{subfigure}{0.2\textwidth}
			\centering
			\begin{tikzpicture}[y=.5cm, x=.5cm,font=\sffamily]

			\draw (0,2) circle (1.35cm);
			\draw (0,2) circle (1.65cm);

			\end{tikzpicture}
			\caption{Shapes}
		\end{subfigure}\label{subfig:shape}
		\begin{subfigure}{0.24\textwidth}
			\centering
			\begin{tikzpicture}[scale=0.5]

			\begin{axis}[
			axis lines = left,
			domain=-1.2:1.2,
			xlabel = x,
			ylabel = y,
			axis equal,
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=-1.2:1.2,
			samples=100,
			color=white,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				-0.61541530  0.656706940 c
				-0.85714932  0.274399412 c
				-0.89116691 -0.125783712 c
				-0.55205903 -0.710795915 c
				0.88053640 -0.186160259 c
				-0.27285906  0.857640913 c
				-0.68979933 -0.578080342 c
				0.89410260  0.102861734 c
				-0.89699390  0.073497901 c
				0.69689897  0.569501381 c
				0.50239031 -0.746728848 c
				-0.80999087  0.392319759 c
				-0.79085117  0.429597991 c
				0.09371138 -0.895107914 c
				0.88025224  0.187499330 c
				-0.44821200  0.780452432 c
				0.38808293  0.812029336 c
				0.89998326 -0.005489085 c
				0.82954861 -0.349068905 c
				0.82433618  0.361206131 c
				-0.76442984  0.475023186 c
				0.85326455  0.286250941 c
				-0.57623985 -0.691337568 c
				-0.86166094  0.259885395 c
				0.31085347  0.844612407 c
				0.46090723 -0.998781520 c
				0.93851411 -0.573751921 c
				-0.09484196  1.095903738 c
				-0.67988031 -0.864732768 c
				-0.59079051  0.927882845 c
				1.07010912 -0.254688989 c
				1.06229667  0.285527214 c
				-0.38382798  1.030861814 c
				0.75791453 -0.797223658 c
				-1.09585907 -0.095356723 c
				0.51360536 -0.972733024 c
				-0.79476398  0.760493405 c
				-0.53630582 -0.960404119 c
				-0.88745589  0.649940029 c
				-0.74717003 -0.807302262 c
				-0.12143215  1.093276833 c
				1.06563952  0.272786398 c
				1.09323023 -0.121850998 c
				0.54446062  0.955804706 c
				0.83251921  0.718965761 c
				0.07220815  1.097627434 c
				-0.36917937 -1.036198143 c
				0.84176190 -0.708122092 c
				-0.59867109  0.922817927 c
				0.99562860 -0.467679040 c
			};
			\end{axis}
			\end{tikzpicture}
			\caption{Scatterplots}
		\end{subfigure}\label{subfig:scatter}
		\begin{subfigure}{0.24\textwidth}
			\centering
			\begin{tikzpicture}[scale=0.5]
			\begin{axis}[
			axis lines = left,
			domain=0:2.1,
			xlabel = Birth,
			ylabel = Death$-$Birth,
			xtick={0,0.4,...,2},
			ytick={0,0.4,...,2},
			scatter/classes={
				c={mark=*,draw=blue}}]
			\addplot [
			domain=0:2.1,
			samples=100,
			color=white,]{x};
			\addplot[scatter,only marks,
			scatter src=explicit symbolic]
			table[meta=label] {
				x     y      label
				0.4950971 1.083311664 c
				0.2349773 0.037147467 c
				0.2256371 0.032775198 c
				0.2311669 0.027387038 c
				0.3926667 0.022625591 c
				0.2545617 0.012075532 c
			};
			\end{axis}
			\end{tikzpicture}
			\caption{PDs}\label{subfig:PD}
		\end{subfigure}
	    \begin{subfigure}{0.29\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/extwo_new}
			\caption{PIs}
			\label{subfig:PI}
     	\end{subfigure}

		\caption{(a) Two shapes, (b) scatterplots of 50 randomly sampled points, (c) transformed one-dimensional persistence diagrams, and (d) persistence images. We use the Gaussian smoothing function with the soft arctangent weight for the persistence images.}
		\label{fig:sim}
	\end{figure}

	We conduct naive two-sample pooled variance $t$-tests at every element to test the differences between the two spaces.
	Figure~\ref{fig:sim.res} presents the average difference of the persistence images, standard errors, and two-sample $t$-test statistics.

	\begin{figure}[!ht]

		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/exdiff}
			\caption{Average difference}
			\label{subfig:ttestdiff}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/exse}
			\caption{Standard errors}
			\label{subfig:ttestse}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=0.95\linewidth]{plots/extstat10}
			\caption{$t$-test statistics}
			\label{subfig:ttestt}
		\end{subfigure}
		\caption{Average differences of one-dimensional persistence images, standard errors, and $t$-test statistics. The $t$-test statistics less than -10 and greater than 10 are plotted as -10 and 10, respectively.}
		\label{fig:sim.res}
	\end{figure}

	As shown in the Figure, some elements which have small mean differences in Figure~\ref{subfig:ttestdiff} correspond to large $t$-test statistics in Figure~\ref{subfig:ttestt} (i.e., elements on the top).
	This happens because the corresponding elements have small standard errors.
	For example, the elements on the top-left corner in Figure~\ref{subfig:ttestt} have large test statistics because the sample standard deviations are relatively small.
	As a result, it is difficult to differentiate whether the large test statistics are due to the mean differences or the small sample standard errors.

	The elements with small mean differences having large test statistics can be seen as the false positives in this case.
	The actual differences of the shapes will appear around the three peak areas in Figure~\ref{subfig:ttestdiff} and \ref{subfig:ttestse}.
	Even though the mean difference of one element is negligibly small compared to other elements, the $t$-test statistic can still be large if its standard error is small enough.

	These false positives are closely related to the vectorization process;
	once persistence diagrams are vectorized, the local element loses the information about the overall variability.
	Therefore, overall mean or variability information should be considered for the hypothesis testing of persistence images.

	Also, when we test multiple hypotheses simultaneously, the multiple testing problem arises. Failure to adjust for multiple testing can lead to an increase of false positives. The multiple testing adjustments such as the family-wise error rate \citep{Hochberg1987} and false discovery rate \citep{Benjamini1995} aim to control the rate of false positives while maintaining statistical power.

	The adjusted multiple hypothesis results need to be robust under dependency among elements in persistence images. Most multiple testing adjustment methods assume independence between tests \citep{Benjamini1995,Storey2002}, whereas some methods control the error rates under various types of dependence \citep{Benjamini2001,Fan2012}. It is necessary to consider potential dependence structures of persistence images in the multiple testing adjustments.

	In the next section, we suggest a hypothesis testing approach for persistent images to handle false positives and adjust multiple tests under dependence.

	\section{Two-stage Hypothesis Test Procedures}
	\label{sec:ts}
	 The suggested hypothesis test approach consists of two stages: filtering and testing.
	 We use two independent statistics (filter and test statistics) in the two stages.

	\paragraph{Stage I: Filtering} The idea of filtering variables has been proposed to increase power for high-dimensional data \citep{McClintick2006,Hackstadt2009,Mieth2016}.
	When the variables are filtered, the number of variables being tested is reduced. As a result, filtering could potentially lead to an increase in the number of discoveries after multiple testing corrections.

	Unlike the way in which filtering has been used in the previous literature, the two-stage procedure in this paper is implemented to drop elements that are noninformative or more likely to be false positives. Therefore, the overall power may not increase after filtering in this instance.

	We add one additional step of removing unnecessary elements in the upper-right corner of the transformed persistent diagram prior to filtering. This part of the image corresponds to an empty area in most persistence diagrams and hence contributes nothing to the inference.

	\cite{Bourgon2010} show that false positive rates are not maintained for a two-stage procedure if inappropriate filter statistics are used.
	More specifically, a filter statistic needs to be independent of the test statistic because the null distribution in the second stage is a conditional distribution given the filter statistic. For example, \cite{Bourgon2010} use the overall sample mean and standard deviation, and \cite{Guo2017} use a sum of squared values as filter statistics for $t$-test statistics.

	In this paper, we use the overall sample mean $\bar{X}^i$ and overall sample standard deviation $s^i$ as filter statistics separately for each element $i$.
	For a given threshold $C$, elements whose filter statistics are less than the $C^{th}$ percentile are removed.

	One of the disadantages of the two-stage hypothesis testing is that the choice of the filtering threshold $C$ could be subjective \citep{Du2014}. \cite{Ignatiadis2016} suggest a data-driven greedy independent filtering procedure that chooses the threshold that maximizes the number of discoveries among all possible candidates. However, this approach is not relevant in our setting because our objective with filtering is to avoid global false positives, not to increase detection power.

	We examine the effect of the choice of filtering threshold in Section~\ref{subsubsec:power} and discuss other possible approaches in Section~\ref{sec:discussion}.

	\paragraph{Stage II: Testing}
	In the second stage, we conduct hypothesis tests on the remaining elements and adjust for multiplicity.
	Though we use the pooled-variance $t$-test in our analysis, any other standard approach can be used such as unpooled $t$-test, nonparametric two-sample tests, ANOVA, and nonparametric ANOVA \citep{Bourgon2010}.

	The effects of the dependence of tests on multiple test adjustment methods  have been studied including positive regression \citep{Benjamini2001}, random correlation structure \citep{Kim2008} and block-correlation structure \citep{Stevens2017}.
	We assume that the elements in persistence images have a block-correlation structure: the elements within blocks have the same correlation coefficients so that the correlation between the elements $i \neq j$,
	\begin{equation*}
	\rho_{i,j}=
	\begin{cases}
	\rho, & \text{if } B_i = B_j
	\\
	0, & \text{otherwise}
	\end{cases}
	\end{equation*}
	where $B_k\in \{1,2,...,B\}$ is the block that the element $k$ belongs to and $B$ is the total number of blocks.
	In persistence images, the elements around peaks can be considered as blocks.
	According to \cite{Stevens2017}, the adaptive procedure of \cite{Benjamini2000}, the q-value method of \cite{Storey2002}, and the two-stage procedure of \cite{Benjamini2006} well control the FDR regardless of the effect size and degree of correlation under the block-correlation settings.
	In this paper, we use the q-value method of \cite{Storey2002}. The q-value is an infimum of a positive false discovery rate (pFDR),
	\begin{equation*}
	pFDR=E(V/R | R>0),
	\end{equation*}
	where $V$ is the number of false positives and $R$ is the number of the hypotheses whose test statistics are in the critical region. The q-values are analogous to the p-values in terms of pFDR \citep{Storey2002}.

	Algorithm~\ref{alg:vec} in Supplementary Material summarizes the suggested testing procedure.

	\section{Results}
	\subsection{Point Cloud Simulation}
	\label{subsec:sim}

	\subsubsection{Effect of Filtering Threshold and Weight}
	\label{subsubsec:power}
	We simulate point clouds from the two shapes of Figure~\ref{fig:sim}.
	First, we randomly draw 10 point clouds of 50 points on shape 1 (one circle) and add noise that follows $N(0,\sigma^2)$.
	Then, we apply the same sampling procedure for shape 2 (two different-sized circles).

	These 20 point clouds, ten from each shape, are used to conduct one hypothesis test.
	We repeat this 500 times for a given $\sigma=0.05$, 0.1, 0.15, and 0.2.

	We construct the Rips complex separately for each of the $20*500*4=40,000$ point clouds and compute persistence homology using the TDA package in R \citep{TDA}.
	The one-dimensional persistence diagrams are converted into persistence images of 40 by 40 elements using a Gaussian kernel with a smoothing parameter $h=0.075$.

	We conduct both the naive and two-stage hypothesis tests.

	\paragraph{Effect of Filtering Threshold}

	We use five filtering thresholds $C=0$\%, 20\%, 40\%, 60\% and 80\% to investigate the effect of $C$. Figure~\ref{fig:sim.thres} shows the q-values of one simulation dataset of $\sigma=0.05$ for the five filtering thresholds when overall sample standard deviations are used as filter statistics. Figure~\ref{Sfig:sim.result} in Supplementary Material shows the q-values when we use the overall sample means as filter statistics.

	Filtering could remove the potential false positive and sparse regions.
	In Figure~\ref{fig:sim.thres}, one of the areas that might be related to false positives is between (1,0) to (2,0); the region that has small mean differences but has small standard errors as Figure~\ref{fig:sim.res}.
	The potential false positive region starts to be filtered out when the threshold reaches 20\%. Also, regions that correspond to sparse areas in the persistence diagrams are filtered out.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/thressd0.pdf}
			\caption{$C=0\%$}
			\label{subfig:0}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/thressd20.pdf}
			\caption{$C=20\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/thressd40.pdf}
			\caption{$C=40\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/thressd60.pdf}
			\caption{$C=60\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/thressd80.pdf}
			\caption{$C=80\%$}
		\end{subfigure}
		\caption{q-values of two-stage hypothesis tests with five threshold values. Overall sample standard deviation is used as a filter statistic.
		}
		\label{fig:sim.thres}
	\end{figure}

	The filtering may not lead to an increase in power, as explained in Section~\ref{sec:ts}. If the filtered noninformative elements correspond to large test statistics, then overall power may decrease.
	Tables~\ref{Stb1}, \ref{Stb2}, \ref{Stb3}, and \ref{Stb4} in Supplementary Material show the power of the two-stage two-sample $t$-tests under $\alpha=0.05$ for the q-value test adjustment procedure, two filter statistics (overall sample mean and standard deviation), five filtering thresholds (0\%, 20\%, 40\%, 60\%, and 80\%), five noise levels (0, 0.05, 0.10, 0.15, and 0.20), and four weights (constant, soft arctangent, hard arctangent, and linear).
	In most cases, positive thresholds have larger powers than no filtering ($C=0\%$). However, when too many elements are filtered out by using a large threshold (e.g., $C=80\%$), then power decreases.

	\paragraph{Effect of Weight}

	We compare four weights (constant, soft arctangent, hard arctangent, and linear) used in the vectorization procedure.
	The latter two (linear and hard arctangent) tend to assign higher weights to persistent features (i.e., points that are far from the 45-degree line in the persistence diagrams). On the other hand, the first two (constant and soft arctangent) give similar weights to features regardless of their persistence.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.49\textwidth}
			\centering

			\includegraphics[width=1\linewidth]{plots/weight_qvalue_mean.pdf}
			\caption{Filter statistics: sample mean}
			\label{subfig:resmean}
		\end{subfigure}
		\begin{subfigure}{0.49\textwidth}
			\centering

			\includegraphics[width=1\linewidth]{plots/weight_qvalue_sd.pdf}
			\caption{Filter statistics: sample standard deviation}
			\label{subfig:ressd}
		\end{subfigure}
		\caption{Powers at significance level $\alpha=0.05$ and filtering threshold $C=50$\%.}
		\label{fig:res}
	\end{figure}

	Figure~\ref{fig:res} compares the powers of two-stage vectorized hypothesis tests with filtering threshold $C=50$\% under these four weights to the permutation test of \cite{Robinson2017}. For the permutation test, pairwise distances between persistence diagrams are measured by the Wasserstein distance.

	We see that the two-stage hypothesis tests with the constant and soft arctangent weights achieve higher power than the permutation tests of \cite{Robinson2017} whereas the linear weight shows lower power. Also, the increase in power over the permutation test is better for higher levels of noise. Similar patterns are seen for either choice of filter statistic.

	The results imply that the proposed two-stage test provides flexible options to compare the differences between collections of persistence diagrams.
	The pairwise distances tend to assign small distances for differences of topological features that are close to the 45-degree line, similar to the linear or hard arctangent weights.
	Therefore, computing the pairwise distances in permutation tests is comparable to measuring differences between persistence images with the linear or hard arctangent weights.

	In our simulated datasets, points close to the 45-degree line in the persistence diagrams play an important role in differentiating the two groups.
	As a result, the permutation tests and two-stage tests using hard arctangent and linear weights do not perform well capturing such differences compared to the tests using the constant and soft arctangent weights.

	\subsubsection{Effect of Persistence Image Resolution}

	We next examine the effect of the resolution of the vectorized persistence diagrams on inference.
	A point cloud of 50 points is sampled from each shape of Figure~\ref{fig:sim} and repeated 50 times so that $n_1=n_2=50$.
	For simplicity, Gaussian noise is not added in this simulation.

	The $50*2=100$ one-dimensional persistence diagrams are converted into persistence images at four different resolutions: 20 by 20, 40 by 40, 60 by 60, and 80 by 80. We use the Gaussian kernel with smoothing parameters set to be 1.5 times the size of the pixel ($h=$0.15, 0.075, 0.05, and 0.0375, respectively). Also, the soft arctangent weight is used. We use overall standard deviation as a filter statistic and pooled-variance $t$-statistic as a test statistic for the two-stage hypothesis tests. The filtering threshold is set to be $C=50$\% and q-values are used to adjust multiple test results.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.40\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/resolution20_sd.pdf}
			\caption{Resolution 20 by 20}
			\label{subfig:rrres20}
		\end{subfigure}
		\begin{subfigure}{0.40\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/resolution40_sd.pdf}
			\caption{Resolution 40 by 40}
			\label{subfig:rrres40}
		\end{subfigure}
		\begin{subfigure}{0.40\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/resolution60_sd.pdf}
			\caption{Resolution 60 by 60}
			\label{subfig:rrres60}
		\end{subfigure}
		\begin{subfigure}{0.40\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/resolution80_sd.pdf}
			\caption{Resolution 80 by 80}
			\label{subfig:rrres80}
		\end{subfigure}
		\caption{q-values at four different resolutions.}
		\label{fig:sim.resolution}
	\end{figure}

	Figure~\ref{fig:sim.resolution} illustrates the filtered elements and q-values for the tested elements. We see that the q-values are similar to each other except for the case of 20 by 20 resolution. These results suggest that the two-stage hypothesis tests are robust to the resolution of the vectorized persistence diagrams unless the resolution is too coarse.

	The persistence image of resolution 20 by 20 does not filter out the elements in the middle, unlike the other higher resolution persistence images.

	The variability of regions in the middle is affected by the neighboring elements because a relatively large smoothing parameter is used for 20 by 20 resolution images.

	\subsection{Two-dimensional Binary Image Simulation}
	\label{subsec:2drock}
	We compare two-dimensional random binary images generated by Algorithm 1 of \cite{Obayashi2018}. Three sets of parameters are used to generate images: $(M=180,S=80)$, $(M=190,S=75)$, and $(M=200,S=70)$, where $M$ is the number of seed points and $S$ is the number of dispersion points. For all settings, we set dispersion scales $\sigma_1=\sigma_2=4$ and threshold $t=0.7$. Figure~\ref{fig:sim.rock} shows examples of generated binary images from each parameter set. The simulated images can be considered as pseudo-rock images. Differences between the images are not easily identifiable by visual inspection.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/N180S80.png}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/N190S75.png}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/N200S70.png}
		\end{subfigure}
		\caption{Examples of two-dimensional pseudo-rock images with parameters $(M=180,S=80)$ (left), $(M=190,S=75)$ (center), and $(M=200,S=70)$ (right).}
		\label{fig:sim.rock}
	\end{figure}

	Hypothesis tests are conducted using topological features of the binary images. We consider three scenarios to compare the pseudo-rock images: scenario 1 examines two groups of $(M=180,S=80)$ images, scenario 2 compares $(M=180,S=80)$ and $(M=190,S=75)$ groups, and scenario 3 tests $(M=180,S=80)$ and $(M=200,S=70)$ groups.

	All groups in each scenario include 50 images of size 200 by 200.

	The images are transformed using a Signed Euclidean Distance Transform (SEDT) procedure \citep{Robins2016}.
	Then, cubical cell complexes are constructed according to the SEDT values.
	Persistent homology is computed using the GUDHI library \citep{gudhi}.
	For the permutation test, $N=200$ repetitions are used.
	For two-stage hypothesis tests, we use the persistence images converted using the four weights and overall mean filter statistic and threshold $C=50\%$.

	The hypothesis test results show that they can differentiate the simulated binary images based on their topological features.

	Table~\ref{tb:rocksim} shows minimum q-values from two-stage hypothesis tests and p-values from permutation tests.

	The case of the largest parameter differences, scenario 3, yields the smallest q-values and p-values.

	\begin{table}[!ht]
		\centering
		\begin{tabular}{|ll|c|c|c|c|c|c|}
				\hline
				&        & \multicolumn{2}{c|}{Scenario 1} & \multicolumn{2}{c|}{Scenario 2} & \multicolumn{2}{c|}{Scenario 3} \\ \cline{3-8}
				&        & Dim 0              & Dim 1             & Dim 0                 & Dim 1                 & Dim 0                 & Dim 1                 \\ \hline
				\multicolumn{1}{|l|}{Two-stage test} & Linear & 0.717              & 0.998             & 0.235                 & 0.135                 & 0.034                 & 0.018                 \\ \cline{2-8}
				\multicolumn{1}{|l|}{}               & Hard   & 0.679              & 0.992             & 0.207                 & 0.277                 & 0.007                 & 0.012                 \\ \cline{2-8}
				\multicolumn{1}{|l|}{}               & Soft   & 0.719              & 0.982             & 0.168                 & 0.220                 & 0.004                 & 0.009                 \\ \cline{2-8}
				\multicolumn{1}{|l|}{}               & Constant     & 0.763              & 0.798             & 0.121                 & 0.792                 & 0.003                 & 0.019                 \\ \hline
				\multicolumn{2}{|l|}{Permutation test}        &   0.520  &   0.700        &    0.090                   &   0.475                    &         0.000              &   0.000                    \\ \hline
			\end{tabular}
	\caption{Minimum q-values of two-stage tests and p-values of permutation tests.}
	\label{tb:rocksim}
	\end{table}

	\subsection{Three-dimensional Sand Pack Image Data Analysis}

	We analyze two sand packs: F42 (unground silica, US Silica Company) and LV60 (Levenseat sand, WBB Minerals, UK).
	These sand pack datasets are obtained by nuclear magnetic resonance scans and micro-computer tomography imaging by \cite{Talabi2009}.
	Each sand pack dataset includes two-samples: F42B and F42C and LV60A and LV60C.
	All sand pack images have $300^3$ voxels with resolution 10.002 $\mu m$.
	Table~\ref{tab:rock} summarizes the sand pack data and Figure~\ref{fig:rock} shows their two-dimensional slice images.

	\begin{table}[H]
		\center
		\begin{tabular}{|c|c|c|c|}
			\hline
			Rock type                         & Sample & Porosity (\%) & Grain surface area ($M^2/m^3$)  \\ \hline
			\multirow{2}{*}{Unground silica F42}

			& F42B   & 33.3 & 44,930     \\ \cline{2-4}
			& F42C   & 33.1 & 45,760      \\ \hline
			\multirow{2}{*}{Levenseat sand LV60}
			& LV60A  & 37.7 & 57,670    \\ \cline{2-4}
			& LV60C  & 37.2 & 61,590      \\ \hline
		\end{tabular}
	\caption{Properties of sand pack data \citep{Talabi2009}.}
	\label{tab:rock}
	\end{table}

	Both rock types have similar porosity (a volume of pores divided by a total volume), but different grain surface area.
	In Figure~\ref{fig:rock}, we see that the unground silica sand pack has larger-sized and circular-shaped grains compared to the Levenseat sand pack. However, it is difficult to compare structural or connectivity differences from the properties or the three-dimensional images.

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.37\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/F42A.png}
			\caption{Unground silica F42}
			\label{subfig:f42}
		\end{subfigure}
		\begin{subfigure}{0.37\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/LV60A_box_new.png}
			\caption{Levenseat sand LV60}
			\label{subfig:lv60}
		\end{subfigure}
		\caption{Two-dimensional slice images of two sand packs. Pores and grains are drawn as white and black elements. The red solid and the blue dotted box areas are two of the subregions that show large variabilities within the same rock image.}
		\label{fig:rock}
	\end{figure}

	We compute persistent homology to extract structural and connectivity information of sand pack data.
	First, we take 27 subregion images of size $100^3$ from each sand pack sample.

	The persistence diagrams are converted into persistence images using the soft arctangent weight.

	We conduct hypothesis tests on three sets of data: between unground silica (F42B and F42C), between Levenseat sand (LV60A and LV60C), and between unground silica and Levenseat sand (F42B and LV60A). We compare two test procedures, the permutation and two-stage hypothesis tests.

	The suggested two-stage hypothesis tests show larger structural and connectivity differences between different types of sand packs compared to the tests between the same type of rocks.
	Figure~\ref{fig:rock.result} shows two-stage hypothesis test results using the overall mean as a filter statistic, pooled variance $t$-tests.
	We see that there are elements that have small q-values for most of the dimensions and rock types.
	However, different types of rocks have larger areas of small q-values than the same type of rocks.
	This implies that various types of topological features account for differences between the two types of rocks.

	We report the other test results using the overall standard deviation in Figure~\ref{Sfig:rock.results} in Supplementary Material.

	\begin{figure}[!ht]
		\centering
		{\hspace{0.7in} F42B and F42C \hspace{0.7in} LV60A and LV60C \hspace{0.6in} F42B and LV60A}\\
		{\rotatebox{90}{Dimension zero}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42_mean_qval_dim0}

			\label{subfig:rockF42dim0}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_LV60_mean_qval_dim0}

			\label{subfig:rockLV60dim0}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42LV60_mean_qval_dim0}

			\label{subfig:rockFLVdim0}
		\end{subfigure}
		{\rotatebox{90}{Dimension one}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42_mean_qval_dim1}

			\label{subfig:rockF42dim1}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_LV60_mean_qval_dim1}

			\label{subfig:rockLV60dim1}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42LV60_mean_qval_dim1}

			\label{subfig:rockFLVdim1}
		\end{subfigure}
		{\rotatebox{90}{Dimension two}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42_mean_qval_dim2}

			\label{subfig:rockF42dim2}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_LV60_mean_qval_dim2}

			\label{subfig:rockLV60dim2}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots/rock_F42LV60_mean_qval_dim2}

			\label{subfig:rockFLVdim2}
		\end{subfigure}
		\caption{q-values of two-sample $t$-tests between F42 (first column), between LV60 (second column), and between F42 and LV60 (third column) for zero- (first row), one- (second row), two-dimensional (third row) persistence images.
		}
		\label{fig:rock.result}
	\end{figure}

	Table~\ref{tb:pt} shows the permutation test results. We use $N=500$ permutations for dimension zero and two, and $N=100$ for dimension one due to the computation time of pairwise distances (See later in this Section for a detailed discussion). Also, we use Wasserstein distance as a pairwise distance. For all cases, the computed permutation p-values are zero; every joint loss of shuffled labels $L(G_{\text{shuffled}})$ is larger than the joint loss of the original unshuffled labels $L(G_{\text{unshuffled}})$. However, the differences of joint losses between unshuffled and shuffled labels are larger for the different types of rocks compared to the same type of rocks.

	\begin{table}[!ht]
		\centering
		\begin{tabular}{|c|c|c|c|c|c|c|}
			\hline
			\multirow{2}{*}{} & \multicolumn{2}{c|}{F42B and F42C} & \multicolumn{2}{c|}{LV60A and LV60C} & \multicolumn{2}{c|}{F42A and LV60C} \\ \cline{2-7}
			& Unshuffled          & Shuffled         & Unshuffled           & Shuffled          & Unshuffled          & Shuffled          \\ \hline
			Dimension 0       & 94.808        & 97.273                &  125.124                 &   131.252               &   121.060               &    207.608              \\ \hline
			Dimension 1       &  222.994                &     228.117            &    274.708               &    329.955              &   254.671               &    528.637              \\ \hline
			Dimension 2       &    135.764              &  137.668               &     175.376         &    183.815              &     158.833          & 263.754                 \\ \hline
		\end{tabular}
		\caption{Joint loss values of permutation tests. For shuffled labels, we report the mean of joint loss values.}
		\label{tb:pt}
	\end{table}

	We want to note that the small q-values and permutation p-values for both methods in the same type of rocks might be due to the small subregion size.
	The subregion size $100^3$ is relatively small to represent the overall rock sample structure.
	The sampled subregions may not have similar structures and connectivity, even when they are taken from the same rock sample.
	For example, we observe that subregions from the same rock sample may have large variabilities themselves in Figure~\ref{subfig:lv60}; the subregion of size 100 by 100 in the top-left corner (in the red solid line box) has smaller sized and more sparse grains than the bottom-right corner subregion (in the blue dashed line box).

	Both test approaches draw similar results, but the proposed two-stage hypothesis test approach has two advantages over the permutation test: interpretation and computational cost.

	First, the two-stage hypothesis test enables detailed inference on the test results; it identifies which areas are the most different in the persistence diagrams. We see the areas with low q-values in Figure~\ref{fig:rock.result} that account for the differences. On the other hand, permutation tests do not provide such information. The permutation test only reports a single permutation p-value that explains overall differences.

	Second, the computational cost for the two-stage test is lower than the permutation test.
	Table~\ref{tb:time} compares the running times of the two approaches for each dimension.
	For all computations, Intel Xeon 2.10 GHz CPUs are used.
	We see that the permutation test requires a substantially longer running time.
	Even if we consider the computation time for the vectorization process as part of the two-stage hypothesis test, the permutation test requires higher computation cost.

	\begin{table}[!ht]
		\centering
		\begin{tabular}{|c|c|c|c|}
			\hline
			\multicolumn{1}{|c|}{Method}                                                       & \multicolumn{1}{c|}{Dimension} & \multicolumn{1}{c|}{Number of cores used} & Computation time \\
			\hline
			\multirow{3}{*}{\begin{tabular}[c]{@{}c@{}}Two-stage\\ test\end{tabular}}    & Dimension 0                    & 1                                         &       $<1$ minute           \\ \cline{2-4}
			& Dimension 1                    & 1                                         &       $<1$ minute           \\ \cline{2-4}
			& Dimension 2                    & 1                                         &    $<1$ minute              \\
			\hline
			\multirow{3}{*}{\begin{tabular}[c]{@{}c@{}}Permutation \\ test\end{tabular}} & Dimension 0 ($N=500$)                    & 40                                        &      $\approx15$ hours           \\ \cline{2-4}
			& Dimension 1 ($N=100$)                   & 50                                        &      $\approx148$ hours            \\ \cline{2-4}
			& Dimension 2 ($N=500$)                 & 40                                        &     $\approx14$ hours              \\  \hline
		\end{tabular}
		\caption{Computation time of two hypothesis test methods for sand pack images.}
		\label{tb:time}
	\end{table}

	Here, the computational bottleneck of the permutation test is the pairwise distance computation.

	For the point cloud simulation in Section~\ref{subsec:sim}, the number of features in a one-dimensional persistence diagram is 13 at maximum.
	On the other hand, we have up to 920 features in a one-dimensional persistence diagram for the sand pack analysis.
	Because we use the Hungarian algorithm to compute the Wasserstein distances ($O(n_f^3)$ time complexity), the computational time is approximately 300,000 times longer for the sand pack analysis.

	Note that we conduct the test on a small number of persistence diagrams, 27 persistence diagrams for each group.
	If we run the permutation test on a large number of persistence diagrams, the computational cost will further increase.

	\section{Discussion}
	\label{sec:discussion}
	This paper presents a two-stage hypothesis test for a vectorized persistence diagram that consists of filtering and testing steps.

	Our approach enables a better interpretation by providing specific regions on the persistence diagram that contribute the most to any observed differences.
	Also, the proposed test procedure is less costly from a computational perspective compared to the permutation test of \cite{Robinson2017}.

	On the other hand, the two-stage test for vectorized persistence diagrams could lead to a potential problem of giving a researcher too many degrees of freedom.
	The proposed method provides a flexible way to conduct a hypothesis test for persistent homology by implementing weights and smoothing.
	The results in Figure~\ref{fig:res} indicate that the weights play an important role for a successful result in revealing the differences between persistence diagrams.
	However, the weight selection process itself is still an unanswered question.
	Also, selecting the filter threshold $C$ can be an issue.
	In practice, one may choose the filter threshold heuristically using an empirical distribution of filter statistics.

	In the future, we plan to explore systematical methods to select the weights and parameters.

	\bibliographystyle{asa}
	\bibliography{References}

	\clearpage
	\section*{Supplementary Material}
	\beginsupplement

	\begin{algorithm}[H]
		\caption{Permutation test \citep{Robinson2017}}
		\label{alg:rt}
		\begin{algorithmic}
			\Require {Persistence diagrams, given group label $G_{\text{unshuffled}}$, number of repetitions $N$, joint loss function $L$}
			\Ensure{Permutation p-value $Z$}
			\State {Compute $L(G_{\text{unshuffled}})$}
			\State {Create an empty vector $L$ of size $N$}
			\For {$i = 1 \to N$}
			\State {Randomly generate shuffled group label $G_{\text{shuffled}}^i$}
			\State {$L[i] \leftarrow L(G_{\text{shuffled}}^i)$}
			\EndFor
			\State {$Z \leftarrow$ sum$(L<L(G_{\text{unshuffled}}))/N$}
		\end{algorithmic}
	\end{algorithm}

	\newpage

	\begin{algorithm}[H]
		\caption{Two-stage hypothesis test procedure}
		\label{alg:vec}
		\begin{algorithmic}
			\Require {Array of $n$ vectorized persistence diagrams $V_{\{m^2 \times n\}}$, given label $G_{\text{given}}=\{I,J\}$, threshold $C$}
			\Ensure{P-values $Z$}
			\State {Create vector $v_x,v_y$ of size $m^2$ that correspond to $x$ and $y$ locations of $V$}

			\State {$V \leftarrow V [ v_x \geq v_y , ]$ }
			\State {Create an empty vector $T$ of size $ \frac{m(m+1)}{2}$}
			\For {$i = 1 \to \frac{m(m+1)}{2}$}
			\State {$T[i]\leftarrow$filter statistic of $i$th element of $V$}
			\EndFor
			\State $t_C \leftarrow C^{th}$ percentile of $T$
			\State $V\leftarrow V[V>t_C,]$
			\State {Create an empty vector $Z$ of size $nrow(V)$}
			\For {$j=1 \to nrow(V)$}
			\State $v_1\leftarrow V[j,I]$
			\State $v_2\leftarrow V[j,J]$
			\State {Conduct a hypothesis test using $v_1$ and $v_j$ and store p-value in $Z[j]$}
			\EndFor
			\State {Apply multiple testing adjustment procedure to $Z$}
		\end{algorithmic}
	\end{algorithm}

	\newpage

		\begin{figure}[!ht]
		\centering

		\begin{subfigure}{0.9\textwidth}
			\centering
			\begin{tikzpicture}
			\draw[thick,->] (-8,0) -- (4,0) node[anchor=south east] {Birth};
			\draw[thick,->] (0,-10) -- (0,4) node[anchor=south west] {Death};
			\draw[thick] (-7.5,-7.5) -- (3,3);
			\node at (0.4,-0.4) {0};

			\draw[thick,dashed] (-6,2) circle (1.8cm);
			\filldraw[fill=black,draw=black] (-4.9,2) circle (14pt);
			\filldraw[fill=black,draw=black] (-6,2) circle (14pt);
			\filldraw[fill=black,draw=black] (-7.1,2) circle (14pt);

			\draw[thick,dashed] (-2,2) circle (1.8cm);
			\filldraw[fill=black,draw=black] (-1.3,2) circle (8pt);
			\filldraw[fill=black,draw=black] (-2,2) circle (8pt);
			\filldraw[fill=black,draw=black] (-2.7,2) circle (8pt);

			\draw[thick,dashed] (-2,-5) circle (1.8cm);
			\draw[thick,blue,->] (-2,-1) -- (-2,-3);
			\filldraw[fill=black,draw=black] (-2.5,-4.3) circle (8pt);
			\filldraw[fill=black,draw=black] (-2,-4.3) circle (8pt);
			\filldraw[fill=black,draw=black] (-1.5,-4.3) circle (8pt);

			\filldraw[fill=black,draw=black] (-1,-5.6) circle (8pt);
			\filldraw[fill=black,draw=black] (-2,-5.6) circle (8pt);
			\filldraw[fill=black,draw=black] (-3,-5.6) circle (8pt);
			\filldraw[fill=black,draw=black] (-3,-5.65) rectangle (-1,-5.55);

			\draw[thick,dashed] (-6,-2) circle (1.8cm);
			\filldraw[fill=black,draw=black] (-5.15,-1.3) circle (14pt);
			\filldraw[fill=black,draw=black] (-6,-1.3) circle (14pt);
			\filldraw[fill=black,draw=black] (-6.85,-1.3) circle (14pt);

			\filldraw[fill=black,draw=black] (-4.9,-2.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-6,-2.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-7.1,-2.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-7,-2.65) rectangle (-5,-2.55);

			\draw[thick,dashed] (-6,-9) circle (1.8cm);
			\draw[thick,blue,->] (-6,-5) -- (-6,-7);
			\filldraw[fill=black,draw=black] (-5.3,-8.3) circle (14pt);
			\filldraw[fill=black,draw=black] (-6,-8.3) circle (14pt);
			\filldraw[fill=black,draw=black] (-6.7,-8.3) circle (14pt);

			\filldraw[fill=black,draw=black] (-4.9,-9.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-6,-9.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-7.1,-9.6) circle (14pt);
			\filldraw[fill=black,draw=black] (-7.5,-9.9) rectangle (-5,-9.3);

			\end{tikzpicture}
		\end{subfigure}
		\caption{Examples of topological features of two-dimensional binary rock images (in dashed circles) and where they appear in the dimension 0 persistence diagram. The dimension 0 persistent diagrams reveal the size of grains and connectivity. For example, larger(smaller) grains have larger(smaller) negative birth values and closely(weakly) connected grains have large(small) negative death values. Also, connected(not connected) grains have negative(positive) death values, so they are recorded on the third(second) quadrant of the persistence diagram. Grains and pores are colored in black and white, respectively.}\label{fig:interpret0}
	\end{figure}

	\newpage

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.9\textwidth}
			\centering
			\begin{tikzpicture}
			\draw[thick,->] (-4,0) -- (8,0) node[anchor=north west] {Birth};
			\draw[thick,->] (0,-2) -- (0,9) node[anchor=south east] {Death};
			\draw[thick] (-1.5,-1.5) -- (8,8);
			\node at (0.4,-0.4) {0};

			\draw[thick,dashed] (-2,6) circle (1.8cm);
			\draw[fill=black,even odd rule]  (-2,6) circle (1.1cm)
			(-2,6) circle (0.6cm);

			\draw[thick,dashed] (2,6) circle (1.8cm);
			\draw[fill=black,even odd rule]  (2,6) circle (1.1cm)
			(2,6) circle (0.6cm);
			\draw[white, fill=white, thick] (2,6.3) -- (2.5,7.5) -- (1.5,7.5) -- cycle;

			\draw[thick,dashed] (9,6) circle (1.8cm);
			\draw[thick,blue,->] (5,6) -- (7,6);
			\draw[fill=black,even odd rule]  (9,6) circle (1.1cm)
			(9,6) circle (0.6cm);
			\draw[white, fill=white, thick] (9,6) -- (10,7.5) -- (8,7.5) -- cycle;

			\draw[thick,dashed] (-2,2) circle (1.8cm);
			\draw[fill=black,even odd rule]  (-2,2) circle (0.7cm)
			(-2,2) circle (0.2cm);

			\draw[thick,dashed] (5,2) circle (1.8cm);
			\draw[thick,blue,->] (1,2) -- (3,2);
			\draw[fill=black,even odd rule]  (5,2) circle (0.7cm)
			(5,2) circle (0.2cm);
			\draw[white, fill=white, thick] (5,2.1) -- (4.6,2.8) -- (5.4,2.8) -- cycle;

			\end{tikzpicture}
		\end{subfigure}
		\caption{Examples of topological features of two-dimensional binary rock images (in dashed circles) and where they appear in the dimension 1 persistence diagram. The dimension 1 persistent diagrams show the structure of grains and pores. The second quadrant (negative birth and positive death) is related to pores trapped inside of grains and the birth(death) value is proportional to the size of the corresponding grain(pore).
		On the other hand, the first quadrant (positive birth and death) represents the grains look like a horseshoe or broken ring. Grains and pores are colored in black and white, respectively.}\label{fig:interpret1}
	\end{figure}

	\newpage

	\begin{figure}[!ht]
		\centering
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/thresmean0.pdf}
			\caption{$C=0\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/thresmean20.pdf}
			\caption{$C=20\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/thresmean40.pdf}
			\caption{$C=40\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/thresmean60.pdf}
			\caption{$C=60\%$}
		\end{subfigure}
		\begin{subfigure}{0.32\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/thresmean80.pdf}
			\caption{$C=80\%$}
		\end{subfigure}
		\caption{q-values of two sample t-tests. The overall sample mean is used as the filter statistic.
		}
		\label{Sfig:sim.result}
	\end{figure}

	\newpage

	\begin{figure}[!ht]
		\centering
		{\hspace{0.7in} F42B and F42C \hspace{0.7in} LV60A and LV60C \hspace{0.6in} F42B and LV60A}\\
		{\rotatebox{90}{Dimension zero}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42_sd_qval_dim0}

			\label{subfig:rockF42dim0s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_LV60_sd_qval_dim0}

			\label{subfig:rockLV60dim0s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42LV60_sd_qval_dim0}

			\label{subfig:rockFLVdim0s}
		\end{subfigure}
		{\rotatebox{90}{Dimension one}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42_sd_qval_dim1}

			\label{subfig:rockF42dim1s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_LV60_sd_qval_dim1}

			\label{subfig:rockLV60dim1s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42LV60_sd_qval_dim1}

			\label{subfig:rockFLVdim1s}
		\end{subfigure}
		{\rotatebox{90}{Dimension two}}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42_sd_qval_dim2}

			\label{subfig:rockF42dim2s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_LV60_sd_qval_dim2}

			\label{subfig:rockLV60dim2s}
		\end{subfigure}
		\begin{subfigure}{0.315\textwidth}
			\centering
			\includegraphics[width=1\linewidth]{plots_supp/rock_F42LV60_sd_qval_dim2}

			\label{subfig:rockFLVdim2s}
		\end{subfigure}
		\caption{q-values of two-sample t-tests between F42 (first column), between LV60 (second column), and between F42 and LV60 (third column) for zero- (first row), one- (second row), two-dimensional (third row) persistence images. We filter using overall standard deviation.
		}
		\label{Sfig:rock.results}
	\end{figure}

	\newpage

	\begin{table}[H]
		\centering
		\begin{tabular}{|c|c|c|c|c|c|}
			\hline
			\multirow{2}{*}{Filter Statistics} & \multirow{2}{*}{Threshold} & \multicolumn{4}{c|}{Noise}            \\ \cline{3-6}
			&       & 0.05  & 0.10  & 0.15  & 0.20  \\ \hline
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Mean\end{tabular}}          & 0\%                          & 0.996 & 0.750 & 0.372 & 0.226 \\ \cline{2-6}
			& 20\%                       & 0.998 & 0.770 & 0.384 & 0.236 \\ \cline{2-6}
			& 40\%                       & 0.998 & 0.796 & 0.432 & 0.266 \\ \cline{2-6}
			& 60\%                       & 0.998 & 0.840 & 0.514 & 0.326 \\ \cline{2-6}
			& 80\%                       & 0.992 & 0.780 & 0.332 & 0.158 \\ \cline{1-6}
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Standard\\Deviation\end{tabular}}       & 0\%                          & 0.996 & 0.750 & 0.372 & 0.226 \\ \cline{2-6}
			& 20\%                       & 0.998 & 0.770 & 0.390 & 0.238 \\ \cline{2-6}
			& 40\%                       & 0.996 & 0.792 & 0.438 & 0.280 \\ \cline{2-6}
			& 60\%                       & 0.994 & 0.728 & 0.450 & 0.310 \\ \cline{2-6}
			& 80\%                       & 0.974 & 0.632 & 0.288 & 0.170 \\ \hline
		\end{tabular}
		\caption{Powers of multiple t-tests on persistence images with constant weight at significance level $\alpha=0.05$}\label{Stb1}
	\end{table}

	\begin{table}[H]
		\centering
		\begin{tabular}{|c|c|c|c|c|c|}
			\hline
			\multirow{2}{*}{Filter Statistics} & \multirow{2}{*}{Threshold} & \multicolumn{4}{c|}{Noise}            \\ \cline{3-6}
			&       & 0.05  & 0.10  & 0.15  & 0.20  \\ \hline
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Mean\end{tabular}}      & 0\%                          & 0.996 & 0.698 & 0.316 & 0.164 \\ \cline{2-6}
			& 20\%                       & 0.998 & 0.712 & 0.322 & 0.168 \\ \cline{2-6}
			& 40\%                       & 0.998 & 0.732 & 0.364 & 0.184 \\ \cline{2-6}
			& 60\%                       & 0.984 & 0.698 & 0.382 & 0.234 \\ \cline{2-6}
			& 80\%                       & 0.856 & 0.486 & 0.352 & 0.226 \\ \cline{1-6}
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Standard\\Deviation\end{tabular}}  & 0\%                          & 0.996 & 0.698 & 0.316 & 0.164 \\ \cline{2-6}
			& 20\%                       & 0.998 & 0.708 & 0.326 & 0.170 \\ \cline{2-6}
			& 40\%                       & 0.992 & 0.712 & 0.374 & 0.212 \\ \cline{2-6}
			& 60\%                       & 0.976 & 0.678 & 0.404 & 0.256 \\ \cline{2-6}
			& 80\%                       & 0.838 & 0.538 & 0.308 & 0.166 \\ \hline
		\end{tabular}
		\caption{Powers of multiple t-tests on persistence images with soft arctangent weight at significance level $\alpha=0.05$}\label{Stb2}
	\end{table}

	\begin{table}[H]
		\centering
		\begin{tabular}{|c|c|c|c|c|c|}
			\hline
			\multirow{2}{*}{Filter Statistics} & \multirow{2}{*}{Threshold} & \multicolumn{4}{c|}{Noise}            \\ \cline{3-6}
			&       & 0.05  & 0.10  & 0.15  & 0.20  \\ \hline
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Mean\end{tabular}}      & 0\%                          & 0.994 & 0.632 & 0.284 & 0.126 \\ \cline{2-6}
			& 20\%                       & 0.996 & 0.646 & 0.296 & 0.136 \\ \cline{2-6}
			& 40\%                       & 0.990 & 0.604 & 0.300 & 0.146 \\ \cline{2-6}
			& 60\%                       & 0.920 & 0.500 & 0.272 & 0.156 \\ \cline{2-6}
			& 80\%                       & 0.824 & 0.414 & 0.220 & 0.124 \\ \cline{1-6}
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Standard\\Deviation\end{tabular}}  & 0\%                          & 0.994 & 0.632 & 0.284 & 0.126 \\ \cline{2-6}
			& 20\%                       & 0.994 & 0.636 & 0.296 & 0.142 \\ \cline{2-6}
			& 40\%                       & 0.944 & 0.554 & 0.300 & 0.164 \\ \cline{2-6}
			& 60\%                       & 0.796 & 0.444 & 0.296 & 0.192 \\ \cline{2-6}
			& 80\%                       & 0.588 & 0.336 & 0.204 & 0.118 \\ \hline
		\end{tabular}
		\caption{Powers of multiple t-tests on persistence images with hard arctangent weight at significance level $\alpha=0.05$}\label{Stb3}
	\end{table}

	\begin{table}[H]
		\centering
		\begin{tabular}{|c|c|c|c|c|c|}
			\hline
			\multirow{2}{*}{Filter Statistics} & \multirow{2}{*}{Threshold} & \multicolumn{4}{c|}{Noise}            \\ \cline{3-6}
			&       & 0.05  & 0.10  & 0.15  & 0.20  \\ \hline
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Mean\end{tabular}}      & 0\%                          & 0.502 & 0.100 & 0.054 & 0.040 \\ \cline{2-6}
			& 20\%                       & 0.356 & 0.124 & 0.068 & 0.042 \\ \cline{2-6}
			& 40\%                       & 0.384 & 0.142 & 0.080 & 0.044 \\ \cline{2-6}
			& 60\%                       & 0.448 & 0.146 & 0.094 & 0.066 \\ \cline{2-6}
			& 80\%                       & 0.532 & 0.132 & 0.092 & 0.056 \\ \cline{1-6}
			\multirow{5}{*}{\begin{tabular}[c]{@{}c@{}}Overall\\Sample\\ Standard\\Deviation\end{tabular}}  & 0\%                          & 0.502 & 0.100 & 0.054 & 0.040 \\ \cline{2-6}
			& 20\%                       & 0.548 & 0.112 & 0.060 & 0.050 \\ \cline{2-6}
			& 40\%                       & 0.564 & 0.118 & 0.062 & 0.058 \\ \cline{2-6}
			& 60\%                       & 0.586 & 0.122 & 0.084 & 0.064 \\ \cline{2-6}
			& 80\%                       & 0.580 & 0.126 & 0.070 & 0.066 \\ \hline
		\end{tabular}
		\caption{Powers of multiple t-tests on persistence images with linear weight at significance level $\alpha=0.05$}\label{Stb4}
	\end{table}

\end{document}


TITLE: The first chern class of Fano manifold
QUESTION [2 upvotes]: If $M$ is a Fano manifold, $L$ is an ample line bundle over $M$. My question is that whether $c_1(L)=\alpha c_1(M)$ for some real number $\alpha$ always holds.

REPLY [3 votes]: No, this is not true. The two classes $c_1(L)$ and $c_1(M)$ live in the abelian group $\operatorname{Pic}(M)$, but for many $M$ this group has rank larger than 1, and there are linearly independent ample classes in the group.
For an example, let $M$ be the blowup of $\mathbf P^2$ in a point. Then $\operatorname{Pic}(M)$ is isomorphic to $\mathbf Z^2$, generated by the hyperplane class $H$ and the exceptional divisor $E$. In this basis $c_1(M)=3H-E$, an ample class. But then for example if $L$ is a line bundle with class $c_1(L)=4H-E$, then $L$ is ample, but the Chern classes are linearly independent.
In fact one can show that ampleness is an open condition on line bundles, so whenever $\operatorname{Pic}(X)$ has rank bigger than 1, we can come up with a similar example.