Topological Data Analysis of Clostridioides difficile Infection and Fecal Microbiota Transplantation
Pavel Petrov Stephen T Rush Zhichun Zhai Christine H Lee Peter T Kim Giseon Heo
q-bio.QM stat.AP 62-07


Topological Data Analysis of Clostridioides difficile Infection and Fecal Microbiota Transplantation
Topological Data Analysis of CDI and FMT 

P Petrov ST Rush Z Zhai CH Lee PT Kim G Heo
Pavel Petrov, MSc Consultant at MNP, 21 West Hastings Street
Vancouver , BC. V6E 0C3 Canada.
Stephen Rush, PhDDepartment of Medical Sciences, University of Örebro, Örebro, 702 81, Sweden. stephen.rush@oru.se
Zhichun Zhai,  PhD Department of Mathematical and Statistical Sciences, University of Alberta, T6G 2G1 Canada.
zhichun1@ualberta.ca
Christine Lee, MD FRCPC. Department of Microbiology, Royal Jubilee Hosptial, Victoria, BC, V8R 1J8, Canada. christine.lee@viha.ca
Peter Kim, PhDDepartment of Mathematics & Statistics, University of Guelph, Guelph, Ontario, N1G 2W1 Canada.  pkim@uoguelph.ca
Giseon Heo, PhDSchool of Dentistry; Department of Mathematical and Statistical Sciences, University of Alberta, Edmonton, T6G 1C9  Canada.
gheo@ualberta.ca




*
    Pavel Petrov  Stephen T. Rush  Zhichun Zhai  Christine H. Lee   Peter T. Kim  and   Giseon Heo
    
==================================================================================================




Computational topologists recently  developed a method, called persistent homology to
analyze data presented in terms of similarity or dissimilarity.  Indeed,
persistent homology studies the evolution of topological features in terms of a single index, and is
able to capture higher order features beyond the usual clustering techniques.    
There are three descriptive statistics of persistent homology, namely barcode, persistence diagram and more recently, persistence landscape.
Persistence landscape  is useful for statistical inference as it belongs to a space  of p-integrable  functions, a separable Banach space.
We apply tools in  both computational topology and statistics to  DNA sequences taken from Clostridioides difficile infected patients treated with
an experimental fecal microbiota transplantation.
Our statistical and topological data analysis are able to  detect  interesting  patterns  among patients and donors.
It also provides visualization of DNA sequences in the form of clusters and loops.



 




Key words and phrases:
Persistent homology; Persistence landscape, Discriminant analysis; Clostridioides difficile; Fecal Microbiota Transplantation; DNA sequences; 16S rRNA gene.










§ INTRODUCTION

Topological data analysis has become a formidable technique for analyzing high-dimensional data, especially when the
purpose is for classification and discrimination.  Methodological advancement has been rampant along with 
applications to medical or scientific data, see for example <cit.>, <cit.> and more recently <cit.>.
In this paper we propose making use of computational topological techniques to analyze gut microbiome data at the basic
DNA sequence level based on data collected from sequencing the 16S rRNA gene.  We are particularly interested in seeing
changes in patient gut microbiome for a certain hypervirulent infectious disease following a radical experimental procedure
that is gaining widespread attention and usage in medicine. 

Clostridioides (formerly Clostridium) difficile (C. difficile) infection (CDI) is the most frequent cause of healthcare-associated 
infections and its rates are growing in the community <cit.>.  
One of the major risk factors for developing CDI is through antibiotics. The 
healthy and diverse bacteria which reside within the colon are the major defense against the growth of C. difficile.  
Antibiotics kill these bacteria and allow C. difficile to multiply, produce 
toxins and cause disease. The current standard of care for this infection are the antibiotics: 
metronidazole, vancomycin and more recently, fidaxomicin.  The efficacy of these antibiotics is limited 
as vancomycin and metronidazole also suppress the growth of anaerobic bacteria such as Bacteriodes 
fragilis group which protect against proliferation of C. difficile.  The efficacy of the recent narrower spectrum fidaxomicin is
still under investigation although the initial data shows promise, <cit.>.
The persistent disruption of healthy 
colonic flora may in part explain the reason for recurrences following a course of treatment with these 
antibiotics.     

As an alternative to antibiotic 
therapy for CDI, in particular for recurrent and refractory diseases, is to infuse healthy gut bacteria directly into the 
colon of infected patients to combat C. difficile by a procedure known as fecal microbiota transplantation (FMT).  
FMT is a process in which a healthy 
donor's 
stool is infused into an affected patient.  This can be performed using
a colonoscope, nasogastric tube, enema, or more recently, in capsulized pill form, <cit.>.  FMT serves to reconstitute the altered colonic flora, in contrast to treatment with antibiotic(s), 
which can further disrupt the establishment of key microbes essential in preventing recurrent CDI. The 
literature reveals a cumulative clinical success rate of over 90% in confirmed recurrent CDI cases <cit.>.

There has been a growing interest into the microbiome of CDI patients <cit.>  especially those 
involved with FMTs <cit.>.  In case of the latter, there are 
differences in:
the route of administration with all forms covered; different donor selection criteria, some used family members, some used a pool of donors; different
sample sizes, although all studies had small sample sizes; and different sequencing procedures and equipment.  Despite these differences, there seems to be two fundamental 
points of agreement across all studies.  The first is that CDI patients have low diversity in their microbiome, and that after receiving an FMT(s), their diversity was increased. The second fundamental agreement is CDI patients who were treated with FMT undergo changes in their microbiome that 
at least initially have similarities to that of their donors.  
This paper provides further reinforcing evidence to support these two fundamental points in the framework
of FMT delivered by enema, with a small exclusive donor pool.  The novelty comes from using computational
topological techniques to demonstrate this.

We now summarize the paper.  In Section <ref> we explain the details of the clinical data.
In Section <ref> we provide topological preliminaries where we go over
the Vietoris-Rips complex, and three topological descriptors: barcodes, persistence diagrams, and persistence landscapes.
In the following Section <ref>, we apply these techniques and demonstrate
the added value of the topological approach at the basic DNA sequence level.  
We complete our article in Section <ref> with a summary of  the key findings. 



§ CDI FMT AND 16S RRNA DATA
 

An earlier report provided details with regard to donor screening and the protocol used <cit.>.
Including the earlier study, FMT was administered to 94 patients and several patients required multiple FMTs to achieve resolution <cit.>.  
From this patient pool we selected 19 patients, not necessarily randomly, for sequencing their 16S rRNA gene prior to treatment, pre-FMT, followed by a
post treatment, post-FMT.  
Consequently, the data that is about to be presented
should not be interpreted with respect to the efficacy of FMT, as the purpose is to better understand the microbial changes that came about as a result of
an FMT(s).

The gender of the patients was 63% female.  The average age (to the time of their first FMT) was 77.11 with a standard deviation of 9.54 years,
range 49 to 92 years. 
In-hospital patients accounted for 53% while the total peripheral white blood count (× 10^9/L) had a mean and standard deviation of 14.6 and 10.87, respectively.  Three patients had temperature greater than 38^o C, and 11 patients experienced abdominal pain.  Approximately half (53%) of the patients were on proton-pump inhibitors.  
Eight patients were refractory to treatment with metronidazole and four patients were refractory to
vancomycin. These four patients were also refractory to metronidazole.  No patients were on fidaxomicin as it had not been available to
the public at that time.
The patient characteristics are provided in Table <ref>, below.




Predisposing conditions that may have resulted in CDI were:  cellulitis, extreme fatigue, respiratory tract infections, septicemia, surgery and open wounds,  and urinary tract infections.  Some conditions were unknown. 
The majority of patients received the following antibiotics prior to contracting CDI: amoxicillin, azithromycin, cefazolin, cefprozil, cephalexin, ciprofloxacin, clindamycin, clarithromycin, cloxacillin, levofloxacin, moxifloxacin and nitrofurantoin.
Some patients claimed no prior antibiotics used.  In an interesting paper, the affects of  ciprofloxacin was studied <cit.> using 16S rRNA deep sequencing. 


The patients all had varying pre-FMT regimens to treat their CDI. All had undergone multiple rounds of traditional antibiotic therapy with metronidazole and/or vancomycin before being administered an FMT. Seventeen received at least one course of metronidazole monotherapy, 18 received vancomycin monotherapy, 6 received vancomycin taper, and 3 received concomitant metronidazole and vancomycin therapy. Patients generally received two courses of metronidazole followed by multiple courses of vancomycin before receiving FMT(s). 
The number of days every patient received each therapy are summarized by their means and standard deviations
and reported below in Table <ref>.






Each patient had two stool samples sequenced, one representing a pre-FMT sequence, and one representing a post-FMT
sequence.
We attempted to sequence each patient prior to them receiving any FMT which for the most part occurred except for one
patient who failed their first FMT and their stool sample sequenced was taken right after that event but prior to their next FMT.
This patient had multiple FMT failures, hence we took that stool sample as the pre-FMT sequence.  We also took a stool sample 
following each patient's pre-FMT sample with at least one FMT in between.  For the most part, the latter occurred following their
last FMT which was the case if the patient resolved their CDI, but if they did not, then the post-FMT sequence was not
necessarily their last.  We would also like to make clear that 4+-FMT means that a patient had at least four treatments but could
have had more.  We indicate such as `4+' to be consistent with our clinical paper, <cit.>.

The breakdown of the data is presented in Table <ref> with the above qualifications.  All patients received
a single
treatment, 1-FMT and 9 of them clinically resolved their CDI.  All who failed the first treatment went on
to receive a second treatment.  Of the remaining patients who received 2-FMT, 2 resolved, while the
remaining went onto to receive a third treatment.  There were 2 successes, and 1 failure, meaning they
did not go on to receive additional treatments.
Of the remaining patients who went on to 4+-FMT, 1 resolved, 1 resolved with antibiotics used in between
treatments, and there were 3 failures.




There were several patients who received antibiotics in between FMTs as described in our previous report <cit.>.
In addition four healthy volunteers served as donors and were screened for transmissible pathogens as was outlined 
in an earlier report <cit.>. The donors took no antibiotics for 6 months prior to stool donation.  Seven donor samples 
taken at various times were sequenced. 

All C. difficile infections were confirmed by in-hospital real-time polymerase chain reaction (PCR) testing for the toxin B gene.
This study sequenced the forward V3-V5 region of the 16S rRNA gene from 19 CDI patients who were treated with FMT(s).  A pre-FMT, a corresponding
post-FMT, and 7 samples from four donors, corresponding altogether to 45 fecal samples were sequenced.
All sequencing was performed on the 454 Life Sciences, GS Junior Titanium Series.
The Qiagen Stool Extraction Kit  (Omega BIO-TEK, Norcross, Georgia) was used to extract the DNA from the fecal samples following the `stool DNA protocol for pathogen detection'.
Subsequent DNA amplification was done using PCR forward and reverse primers.
One round of DNA amplicon purification was performed using the QIAquick PCR Purification Kit (Qiagen, Valencia, CA) followed by two rounds of purification using Agencourt AMPure XP beads (Beckman Coulter Inc., Mississauga, ON).  

We examined
19 pairs of pre-FMT and post-FMT patients, as well as donors, selected from the 94 CDI patients treated by the first author over the period 2008-12, <cit.>. 
The selection was not random and was chosen to reflect as wide a variability as possible.  Furthermore, we deliberately chose 4 FMT failures to try and
better understand what FMT can and cannot do.  It was perhaps through the failures that we learned the most.  The failed FMT cases were ultimately
resolved with antibiotics even though some patients were refractory to metronidazle and/or vancomycin beforehand.  
Indeed FMT acted as a `gut primer' for antibiotics to fulfill it's role in clearing recurrent and refractory CDI. 
Although this phenomenon of
failed FMT patients resolved with antibiotics afterwards has been observed in our work, <cit.>, as well as others, see for example <cit.>, as far as we are
aware this is the first paper that has pursued this from a microbiome point of view.
The data comes from the clinical work of the fourth author who is a practising physician.
Some of this data
was also examined in <cit.> and <cit.>.

Table <ref> below, provides some descriptive statistics about the total and unique number of DNA sequences found in the 45 samples. 
Note that it is impossible to ensure that we have an approximately equal number of sequences in each sample <cit.>. 



Consider the DNA sequences λ and μ. Let x be the number of point dissimilarities between DNA bases in λ and μ, and let x_λ and x_μ be the number of DNA bases in λ and μ, respectively. Let x_O be the number of places where one sequence contains a DNA base and the other an O, a gap. Let y_λ and y_μ be the number of gaps in λ and μ, respectively. Let z be the length of the alignment sequences; we can assume they are
equal otherwise we choose the smaller.  The various dissimilarity metrics are the following.

The one-gap distance d_O is defined by 

    d_O(λ,μ)=x+min{y_λ,y_μ}/min{x_λ+y_λ,x_μ+y_μ}.

It treats a string of O's flanked by any DNA bases in one sequence and the corresponding region in the other as one mismatch.

The no-gap distance d_N is defined by 

    d_N(λ,μ)=x/min{x_λ,x_μ}.

It ignores gaps in a sequence and the corresponding region in the other sequence.

The each-gap distance d_E is defined by 

    d_E(λ,μ)=x+x_O/z.

It treats every DNA base-placeholder pair as a mismatch.

These metrics produce true mathematical distance matrices, which are symmetric and contain zeros on the diagonal.
For this and most studies, the metric of interest is the one-gap metric, <cit.>.  As an example, suppose there are two DNA sequences with the following base pair orientation:
 and  .
Here there are two mismatches and one gap. The distance is calculated as the number of mismatches divided by length of the shorter sequence. The length of the shorter sequence is 10 base pairs, since the gap is considered a single position. Hence distance is 0.3, <cit.>.
We can use the other definitions of distance, but the definition provided is the one most commonly used. In addition, the results do not vary significantly. 
Most of the DNA sequences found in the samples appear several times, and hence the distance between them will be zero as they are identical. For this reason, the unique DNA sequences are taken.

Table <ref> shows the minimum number of unique sequences is 147. DNA sequencing does not provide exact results, hence as the number of sequences in a sample increases, some mutations inevitably occur and these are recorded as unique sequences <cit.>. In other words, as the total number of sequences increases, the number of unique sequences is also inflated. For this reason, it is necessary to subsample from the number of unique sequences. The smallest number of unique sequences is 147, hence a weighted subsample of size 147 is taken from the number of unique sequences for our research. 
The pairwise distance between the 147 sequences is calculated in each of the 45 samples using the one-gap metric.
Thus the data used for the primary analysis is 147 × 147 distance matrices for each of the 45 patients and donors. 








§ SOME TOPOLOGICAL PRELIMINARIES


Persistent homology is a branch of computational topology that has been popularized by <cit.> and 
<cit.>.
Several researchers have shown that persistent  homology  works well on detecting topological and geometrical  features in high dimensional data, see <cit.> and  <cit.>,  for example. 
Let us suppose that we have points on a manifold whose dimension is not necessarily known. 
Topology   studies the connectivity of these points.
Each point  is replaced by a disk (ball) with  radius  centered at each point.  
If the disks overlap, connect  those  center points with edges.  As the radius increases, more points will be connected and so the number of connected components will decrease. This is analogous to clustering analysis in statistics. The connected components (clusters) are considered as 0-degree topological features.
The number of connected components  is denoted as the 0-th Betti number, β_0.
Connecting the points with edges will also create simplices (convex hull of a geometrically independent set of points) and  produce topological features in higher dimensions. The low dimensional simplices are well known;   a vertex (0-degree simplex), an edge (1-degree simplex), a triangle (2-degree simplex), and a tetrahedron (3-degree simplex).
A loop  is also called a 1-degree  topological feature  and a void is a 2-degree topological feature.
The k-th Betti number β_k counts the number of k-degree topological features (k-dimensional `holes').

Persistent homology studies the history of topological features as the parameter ε varies.
It records  the time at which a topological feature appears and disappears.
The birth, death, and  survival time of features are recorded as a barcode <cit.>.
A true feature in the data  lives over a long time while noise is short lived.
We illustrate a fundamental idea of persistent homology  in Figure <ref> with 40 randomly selected points from a double  annulus.  The loops in the  middle of each annulus  are prominent  1-degree features  and their  intervals in the barcode shows their persistence.





Clusters  come in different shapes, see a few synthetic examples in Figure <ref>.
Clusters are  homogeneous subgroups where the meaning of homogeneity is dependent on the types of similarity measure. 
In the double annulus, Figure <ref>, if we consider geodesic distance between points,  those points are connected to form simplices (edges and triangles).  All the edges and triangles together form a band which is homotopy equivalent to a circle (1-degree topological feature).  The authors in <cit.> applied support vector clustering to concentric rings and were able to detect three rings as clusters. In Figure <ref>, we generate similar concentric rings as in <cit.>.  Persistent homology analysis shows three persistent loops which correspond to three clusters in <cit.>. Consider another example,  where points are randomly sampled from  the number 8 with noise, see Figure <ref>. 
What are the clusters in this data?  If one thinks of `blobs' as clusters, there are several patches made up with a few points.
However,  each loop, 1-degree topological feature,  is represented as the set of points which are close in the sense of the shortest path distance.




To study the connectivity  of a space,  the space need not  be  represented  as a point cloud. 
All we need to know is how close  the points are in a space where they live.  
Thus point clouds or matrices of (dis-)similarity measurements among other data forms are input data for analysis of persistent homology.  
A good reference for algebraic  topology and persistent homology is  <cit.>. 









 §.§ Vietoris-Rips complex and topological descriptors 

We will further explain Figure <ref> and then introduce three topological descriptors, barcode, persistence diagram <cit.> and persistence landscape <cit.>.
A collection of point cloud data in a metric space is converted to  a combinatorial graph whose edges are determined by closeness between the points. 
While a  graph captures  connectivity and clustering of data, it  ignores higher dimensional features.
This  idea of graph can be extended to a simplicial complex, which is a collection of simplices.

Suppose there is a finite set of points {v_i}_1^n  in  and (dis-)similarity measure is denoted as d(v_i, v_j).
A k-simplex is  a set of all points x∈ such that x=∑_i=1^k, a_iv_i, where ∑ a_i=1, a_i ≥ 0.
It is easy to picture low dimensional simplices; a 0-simplex a vertex, a 1-simplex an edge joining two vertices,  a 2-simplex a triangle, 
 a 3-simplex a tetrahedron.
Vietoris-Rips complex _ is a set of simplices whose vertices have pairwise  distance within d(v_i, v_j) ≤.
Algebraic topology adds group structure onto the complex,  H_k(_), called k-th homology group. 
For coefficients in a field, H_k(_) is a vector space, whose basis consist of linearly independent k-dimensional cycles that are not  boundaries. 
The k-th Betti number β_k, is the rank of H_k(_) and  counts the number of k-dimensional holes of a simplcial complex. 

At each fixed , homology group H_k(_)  can be calculated, but  computational topologists think of persistence.  The simplical complexes grow as  increases, that is,  𝒱_⊆𝒱_*, for ≤^*, this inclusion induces a linear map,  H_k(_) → H_k(_^*). This allows us to examine the filtration of homology.   
The evolution of the simplicial complexes over  increasing values of  can be completely  tracked using barcodes or  persistence diagrams.  
Barcode is the  multiset of intervals (', ”), where '  and ” indicate birth and death time of a  topological feature. 
Alternatively,  the birth and death times can be represented by a point (',”) in ℝ^2. The collection of these points in ℝ^2 is a persistence diagram (see Figure <ref>). Vietoris-Rips complex and its evolution are demonstrated with random points selected from a double annulus, see Figure <ref> above. 

In a simple summary, the point cloud data are transformed to barcodes or  persistence diagrams in each  dimension.
Is it possible to calculate means and variances of barcodes or persistence diagrams?
It is well known that  the Fréchet mean of barcode (persistence diagram) is not unique.  
Many researchers  have advanced in developing  theories  in this research area <cit.>.  
Bubenik <cit.> introduced  a third topological descriptor, persistence landscape.  
Given an interval (b, d), with b≤ d, define a function leading to an isosceles triangle,  f_(b,d): →, 
f_(b, d)=min (t-b, d-t)_+,  where u_+=max(u,0). 
The persistence landscape  corresponds to a multiset of intervals {(b_i, d_i) : b_i ≤ d_i} and to a set of functions, 
{λ(k, t): ×→} where : λ (k, t) is the  k-th largest value of {f_(b_i, d_i )}.
See an illustration in Figure <ref>.











 §.§  Statistical Inference with Persistence Landscape
 
On ×, the product of the counting measure on  and the Lebesgue measure on  are used.
The persistence landscape,  λ(k,t): ×→,  is bounded and nonzero on a bounded domain.
Hence persistence landscape belongs to ^p(×), with a  metric induced by p-integrable functions and hence  is a separable Banach space <cit.>.
Bubenik also showed that when p ≥ 2, with finite first and second  moments, persistence landscape satisfies a Strong Law of Large Numbers (SLLN) and a Central Limit Theorem (CLT).
Suppose Λ_1, …, Λ_n are the random variables corresponding  to persistence landscapes. The vector space structure of ^p(×) induces  the mean 
landscape as the pointwise mean, λ̅(k,t)=1/n∑_i=1^nλ_i(k, t).

On p-integrable space,  ^p(×), where p ≥ 2  and under the assumption of  finite first and second moments, 
for any continuous linear functional f, the random variable f(λ(k,t)) also satisfies  SLLN and CLT <cit.>.
There are many choices of f, but the integration of λ might be a natural choice,
f(λ(k,t))=∑_k∫_λ(k, t) dt. It has a good interpretation; the  values  of f 
are an enclosed total  area of all curves  λ(k,t). Choice of f might depend on  data, see other choices in <cit.>.   
With the integration as a functional choice, 
we are ready to set hypotheses.
Let Y_1=f(λ_1(k,t)) and Y_2=f(λ_2(k,t)) be random variables for groups 1 and 2.
We let μ_1 and μ_2 be corresponding population means. The hypothesis of interest is

    H_0: μ_1-μ_2=0       H_a: μ_1- μ_2 ≠ 0.




We also consider the similarity measure between persistence landscapes  as well as between persistence diagrams.
The measure between persistence landscape is defined as the p-norm of difference.
Suppose there are two samples that have landscapes denoted as λ_1(k, t) and λ_2(k, t). Then the L_p-norm is defined in (<ref>). This compares pairwise the area under the contours between λ_1(k,t) and λ_2(k,t)




    ||λ_1-λ_2||_p=(∑_k∫_ℝ|λ_1(k,t)-λ_2(k,t)|^p)^1/p.

For our data, we will calculate the persistence landscape  distance using both L_1 and L_2 norm in  (<ref>).
Wasserstein distance (see  (<cit.>, for example) is a popular  measure of dissimilarity  between persistence diagrams.
All the results based on persistance landscape (p=1, 2)  and Wasserstein distance are similar and so all our statistics and presentations are  based on L_2 norm in  the following sections. 







§ TOPOLOGICAL DATA ANALYSIS


The computation of the Vietoris-Rips complex, which is computationally intensive, is  carried out using the  package in R <cit.> on the Westgrid computer network. 
All 45 samples had intervals in barcode in degrees zero and one. However in degree two, only 12 of the pre-FMT samples and 15 of the post-FMT samples have intervals in barcode. Persistence diagram  makes visual pairwise comparisons easier than barcode. 
Figure <ref> shows the corresponding persistence diagrams of patient 10 before and after the FMT  treatment. 



At first glance there does not appear to be anything interesting but there are some general trends. We observe that the degree 0 persistence  diagram for the pre-FMT samples have components that die a lot sooner than those in the post-FMT samples. Similarly, the birth and death times of these loops are shorter for the pre-FMT samples than for the post-FMT samples.
These observations  may indicate that there could be a true difference in the topological structure in degrees zero and one between the two groups. On the other hand,  many of the intervals in degree one  are very short, thus what we observe may be noise rather than signal. 

Unlike barcode or persistence diagrams, we can calculate means and variances of persistence landscapes.
Figure <ref> below, shows the average persistence landscapes of pre-FMT and post-FMT samples in degree 0 and 1. 
The trend in average persistence landscapes is the same as in persistence diagrams or barcodes, that is, the post-FMT samples have slightly longer intervals than the pre-FMT in degree 0 and 1.



Quadratic discriminant analysis (QDA)  was performed on  the β_0- and β_1-Isomap embedded coordinates.
Three groups are well separated on the plane of both degree 0 and 1. 
It is interesting to see that the donors are grouped on one side and pre-FMT patients  on the other side, and  post-FMT patients between the two. 
This is in complete agreement with what clinicians believe is happening and it is frequently reported that patients gut microbiome take on the characteristics of the donor microbiome
following an FMT, see <cit.>. 

It would have been very interesting to compare the matched donor and the patient after FMT, but this information was not recorded during this study, <cit.>.
The classification in ℝ^2 space  in degree 0  in Figure <ref> shows the patients 7, 16 and 19 post-FMT, become much like donors, particularly  donors 2 and 3. Patients  9 and 15 post-FMT  did not seem improved as they are similar to the cases pre-FMT.
The classification in ℝ^2 space  in degree 1  in <ref>  show a few post-FMT patients  (7, 10, 15) are close to the donors after FMT while most post-FMT patients remain similar to  pre-FMT.  Following <cit.>, a clinical trial comparing the efficacy of frozen versus fresh FMT has been completed, see <cit.>.  Here stool samples were
collected at pre-FMT, followed by
day-10, week-5, week-13 following a patients last FMT along with the exact donor stool sample pairing.  Sequencing of this data is currently underway.





We recall that our `raw' data was the dissimilarity matrix between 147  uniques DNA sequences per subject, so a total of forty five 147× 147
matrices. For each matrix, we construct a Vietoris-Rips complex, then calculate  persistence landscapes which enables us to perform statistical inference. 
Hypothesis  test (<ref>) were carried out to comparing pre-FMT and post-FMT samples. 
The p-values of the paired t-tests are 0.0064, 0.0083 and 0.2591 in degree 0, 1, and 2, respectively.
Since all the analysis above is based on 147 sequences randomly chosen from  those patients whose number of  of DNA sequences are bigger than 147,   we repeated analysis 10 times with  independent 147 DNA  samples. The test statistics  and p-values for 10 runs are similar showing consistency of result regardless of which  147 DNA sequences were applied. 

Hypothesis tests  show significant difference  of topological features in degree 0 and 1 between patients before and after FMT treatment.  
For degree 0, this implies that  the number of  clusters  and their persistence on DNA sequences in pre-FMT samples  are different from  those in post-FMT samples. For degree 1, this implies that  the  number  of loops (cycles) and their persistence in DNA sequences in pre-FMT samples  are different from  those post-FMT samples.
We  present  DNA sequences  as points  cloud in ^3 and observe  patterns of clustering and loops in DNA sequences in  the following section.






 §.§ Clusters and loops in  DNA sequences


Applying dimensional reduction methods; Isomap and multidimensional scaling (MDS), to dissimilarity measure between DNA sequences, we project the DNA sequences to ^3 and obtain embedded coordinates.
Scree plots appear to indicate embedding dimension of the DNA sequences is 3.
The  residual variance  in MDS was much higher than for Isomap, hence figures based on Isomap  are presented below.






Figure <ref> shows the Isomap embedded coordinates for donor 3, post-FMT patient 7 and pre-FMT/post-FMT patient 19.
 From the Figure <ref>, it can be seen that donor 3 and post-FMT patient 7 have a similar spread among the sequences and there are 2-3 distinct clusters. 
We have also  noticed that donor 2 and post-FMT patient 18 also have similar structures. 
The donors generally have a wider spread, which would indicate more diverse DNA sequences and hence  a healthier gut microbiome.
For the number of of clusters, there are higher number of clusters in post-FMT samples,  however  some clusters contain only one DNA sequence  and  are more spread out.
For example, there are about 5 groups in  pre-FMT patient 19;   many clusters with fewer DNA sequences in each cluster in post-FMT patient 19.  
Ignoring the singleton clusters, we observe two  large clusters in post-FMT patient 19.
 This information was  shown  in Figure <ref> as well as on  barcode and persistence diagram.  The two clusters in post-FMT patient 19 can be seen in the barcode diagram as the two points that have the highest `death' time. The spread out clusters in pre-FMT 19 are shown as the bars that have the earlier birth and earlier death times
(figures are not shown here). This trend is visible in other samples and this might explain the small p-value for testing the difference in area under the persistence landscapes which was calculated in Section <ref>.











§ CONCLUSIONS
 
We illustrated how topological data analysis can be applied to similarity measures data of DNA sequences.
DNA sequence data was analyzed using three summary statistics of persistent homology; namely, barcodes, persistence diagrams and persistence landscapes. The main objective was to see if there are any differences in the topological feature of DNA sequences in the gut microbiome of CDI patients before and after FMT treatment. 

From visual inspection of barcodes and persistence diagrams it was seen that the components in dimensions zero and one died sooner in the pre-FMT samples than in the post-FMT samples. Persistence landscapes were able to present this difference more formally, showing that there is a difference in the average area under the persistence landscapes for pre-FMT and post samples. Alternative interpretation of this was that there was a difference in the number and size of clusters in dimension zero, and the number and size  of loops in dimension one. The post-FMT  samples had more clusters than the post  samples, whereas there was no visually obvious difference in the size of the loops, but the loops were bigger for the post samples.

We performed discriminant analysis on β_0, β_1 -Isomap embedded coordinates.
The classification on two dimensional space for both degree 0 and 1  show good separation among three groups.
For degree 0, the first two  discriminant function  separates  pre-FMT  patients and donors, and post-FMT patients sit  between the two.
For degree 1, the first two  discriminant function  separates  pre-FMT and post-FMT patients and donors are between  the two.


The major drawbacks of this analysis is that information about individual sequences is lost. This project looked at the topological structure created by the sequences but no details were provided about the individual sequences. As the micro-biological technology and methods improve it may be interesting to incorporate this information.  Also of interest would be meaningful identification of bacterial species. Studies of the 16S rRNA gene only measure presence of the species, but do not say anything about their functionality.  For the latter one has to turn to the metabolome which is currently
under investigation with the patients in <cit.>.



We are grateful to all of the donors, families and patients 
who took part in this study.  We also appreciate the clinical and research staffs at St Joseph's Healthcare Hamilton
where the clinical work had been performed.  The corresponding
author would also like to thank the participants of the SAMSI Working Group “Nonlinear 
Low-dimensional Structures in High-dimensions for Biological Data" which was part of the 2013-14 SAMSI LDHD Program.  Much of the
discussions were centred on the work presented.

We  also thank Violeta Kovacev-Nikolic for her help with matlab code and Figure 1;  Professor Patrick Schloss for his help using mothur; and Yi Zhou for his help with Figures 2 and 3.
Computations in this research were largely enabled by resources provided by WestGrid and Compute Canada.

We would like to acknowledge funding support provided by:  CANSII CRT;  CIHR 413548-2012; McIntyre Memorial Fund; Michael Smith Foundation; NSERC DG 293180, 46204; NSF DMS-1127914; and, PSI Foundation Health Research Grant 2013, 2017.

The study and permission protocol was approved by the Hamilton Integrated
Research Ethics Board #12-3683, the University of Guelph
Research Ethics Board 12AU013 and the University of Alberta,  Health Ethics approval Pro00047221.











abbrvnat





