1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375
|
\documentclass[11pt]{article} % use larger type; default would be 10pt
\usepackage[utf8]{inputenc} % set input encoding (not needed with XeLaTeX)
%%% Examples of Article customizations
% These packages are optional, depending whether you want the features they provide.
% See the LaTeX Companion or other references for full information.
%%% PAGE DIMENSIONS
\usepackage{geometry} % to change the page dimensions
\geometry{a4paper} % or letterpaper (US) or a5paper or....
% \geometry{margin=2in} % for example, change the margins to 2 inches all round
% \geometry{landscape} % set up the page for landscape
% read geometry.pdf for detailed page layout information
\usepackage{graphicx} % support the \includegraphics command and options
% \usepackage[parfill]{parskip} % Activate to begin paragraphs with an empty line rather than an indent
%%% PACKAGES
\usepackage{booktabs} % for much better looking tables
\usepackage{array} % for better arrays (eg matrices) in maths
\usepackage{paralist} % very flexible & customisable lists (eg. enumerate/itemize, etc.)
\usepackage{verbatim} % adds environment for commenting out blocks of text & for better verbatim
\usepackage{subfig} % make it possible to include more than one captioned figure/table in a single float
% These packages are all incorporated in the memoir class to one degree or another...
%%% HEADERS & FOOTERS
\usepackage{fancyhdr} % This should be set AFTER setting up the page geometry
\pagestyle{fancy} % options: empty , plain , fancy
\renewcommand{\headrulewidth}{0pt} % customise the layout...
\lhead{}\chead{}\rhead{}
\lfoot{}\cfoot{\thepage}\rfoot{}
%%% SECTION TITLE APPEARANCE
\usepackage{sectsty}
\allsectionsfont{\sffamily\mdseries\upshape} % (See the fntguide.pdf for font help)
% (This matches ConTeXt defaults)
%%% ToC (table of contents) APPEARANCE
\usepackage[nottoc,notlof,notlot]{tocbibind} % Put the bibliography in the ToC
\usepackage[titles,subfigure]{tocloft} % Alter the style of the Table of Contents
\renewcommand{\cftsecfont}{\rmfamily\mdseries\upshape}
\renewcommand{\cftsecpagefont}{\rmfamily\mdseries\upshape} % No bold!
%\usepackage{url}
\usepackage{xspace}
\usepackage{hyperref}
%%% END Article customizations
%%% commands
\newcommand{\sepp}{SEPP\xspace}
\newcommand{\ins}[1]{{\tt #1}}
\newcommand{\seppfile}{sepp.tar.gz\xspace}
\newcommand{\file}[1]{{\sf #1}}
\newcommand{\sate}{SAT\'{e}\xspace}
\newcommand{\pplacer}{pplacer\xspace}
\newcommand{\guppy}{guppy\xspace}
\newcommand{\hmmer}{HMMER\xspace}
\newcommand{\raxml}{RAxML\xspace}
\newcommand{\arch}{Archaeopteryx\xspace}
\newcommand{\sbun}{\textasciitilde/.sepp/bundled-v3.1}
\newcommand{\sepphome} {\file{{\raise.17ex\hbox{$\scriptstyle\sim$}}/sepp}\xspace}
%%% The "real" document content comes below...
\title{A brief tutorial on \sepp }
\author{Siavash Mirarab, Nam Nguyen, Tandy Warnow}
%\date{} % Activate to display a given date or no date (if empty),
% otherwise the current date is printed
\begin{document}
\maketitle
\tableofcontents
\section{Introduction to \sepp}
\sepp stands for \sate-Enabled Phylogenetic Placement~\cite{sepp} and addresses the problem of phylogenetic placement for meta-genomic short reads. More precisely, \sepp addresses the following problem.
\begin{description}
\item [Input:] {\em backbone} tree T and alignment A for a set of full-length gene sequences, and set X of fragmentary sequences for the same gene
\item [Output:] placement of each fragment in X into the tree T, and alignment of each fragment in X to the alignment A.
\end{description}
Phylogenetic placement puts unknown short fragments into a phylogenetic context, and hence helps identifying species included in a metagenomic dataset.
Phylogenetic placement involves two steps: alignment of short fragments to full length sequence alignment A (also called the {\em backbone} alignment), and then placement of aligned short reads on the fixed tree T (also called the {\em backbone} tree).
\sepp operates by using a divide-and-conquer strategy adopted from \sate \cite{liu09science,sate2} to improve the alignment of fragments to the
backbone alignment (produced by running \hmmer \cite{hmmer}).
It then places each fragment into the user-provided tree using \pplacer \cite{pplacer}.
Our study shows that \sepp provides improved accuracy for quickly evolving genes as compared to other methods. For more information see \cite{sepp} (available at: \url{http://www.ncbi.nlm.nih.gov/pubmed/22174280}).
\section{Installing \sepp}
First you are going to setup \sepp on your machines. \sepp currently runs only on Linux and Mac. Those running a Windows need to either use cygwin (\url{http://www.cygwin.com/}), or a Ubuntu virtual machine image we have created.
\subsection{Installing on a Linux Machine}
To download and install \sepp, follow these steps:
\begin{enumerate}
\item Download software from \url{https://github.com/smirarab/sepp/zipball/master}.
\item Copy the archive to your favorite location (we will use \sepphome in this tutorial), and unpack the zip file using your favorite software.
\item Open a terminal and change into the unpacked directory (\ins{cd \sepphome}).
\item Run the following command:\\
\ins{python setup.py config}\\
Note that this step needs to be run by each individual user. This step creates a \file{.sepp/main.config}
file in your home directory and copies the bundled tool under the same directory.
This file contains defaults settings of sepp, and can be modified by users if needed.
If you wish not to use the home directory, you can use
\ins{python setup.py config -c}\\
instead. This will ensure that the main config file and binary files are written to your current
directory instead of the home.
\item In the new directory, there should be a \file{setup.py} file. Run the following command:\\
\ins{sudo python setup.py install}\\
If you don't have root access, instead use the following command (or use \ins{--prefix} option, as described in the \file{README} file)\\
\ins{python setup.py install --user}
\end{enumerate}
Refer to the README file for more information regarding the installation, and solutions to common installation problems.
\subsection{Using Virtual Machine}
An Ubuntu VM image with \sepp installed on it is available for download at \url{http://www.cs.utexas.edu/~phylo/software/Phylolab.ova}.
If you were unable to install \sepp on your machine, you can use this VM image. You first need to copy the VM image to your machine. Then, open the VM image in your favorite VM software (e.g. VirtualBox) and start the VM. Once your VM starts, you can run \sepp from a terminal. The username and password for the virtual machine are ``ubuntu'' and ``reverse''.
\section{Running \sepp}
\sepp is currently available only as a command line tool, and
so the tutorial is based on this command line usage.
In this section we will run \sepp on a small sample dataset provided with the software.
The sample dataset consists of a \sate backbone alignment and tree on the ``pyrg'' marker gene with only 65 sequences
(previously studied in \cite{metaphyler}).
The fragments are from a WGS sample of a mock community created by the NIH Human Microbiome Project (\url{http://www.hmpdacc.org/HMMC/}).
The fragment file we provide includes only 106 fragments that we found to possibly belong to ``pyrg'' marker (based on hits to its HMMER profile).
\subsection{Run \sepp with -h option to see the help }
\begin{itemize}
\item Make sure you have a terminal open.
\item Make a directory where you will run \sepp and cd into it (e.g. \ins{mkdir seppRuns; cd seppRuns;} )
\item run \sepp with -h option to see a help:\\
\ins{run\_sepp.py -h}
\end{itemize}
\subsection{Sample Datasets - Default parameters}
\begin{itemize}
\item Copy test datasets into your directory. Test datasets are part of the distribution and can be found under \file{\sepphome/test/unittest/data} (use the directory where you unpacked \sepp instead of \sepphome).
For example, on Unix run:\\ \ins{cp -R \sepphome/test/unittest/data/* .}\\
Note that on the VM, the files to copy are found in \textasciitilde/testdata/sepp/.
\item Execute the following command to run \sepp on a sample biological dataset.\\
\ins{run\_sepp.py -t mock/pyrg/sate.tre -r mock/pyrg/sate.tre.RAxML\_info -a mock/pyrg/sate.fasta -f mock/pyrg/pyrg.even.fas} \\
(\sepp should finish running in about 2 minutes.)
\end{itemize}
This runs \sepp on the given example dataset. All the options provided to \sepp in this example run are mandatory. As you can see, there are few inputs required to run \sepp. The following describes the minimum input of \sepp:
\begin {description}
\item[Backbone tree:] this is the tree on which \sepp places short fragments.
This tree should be a binary maximum likelihood (ML) tree
in newick format; we therefore recommend you estimate
the ML tree using \raxml~\cite{raxml} or phyml~\cite{phyml}\footnote{Note - we
haven't tested \sepp with phyml trees yet, and all our analyses are
based on RAxML.} on the backbone alignment (see below). The
input tree is given to \sepp using -t option.
\item[Backbone Alignment:] this is a multiple sequence alignment
of full length sequences for some gene. These sequences need to
be for the same gene as the fragmentary sequences, as phylogenetic
placemeent only makes sense in this case. The backbone
alignment needs to be highly accurate, since it determines the backbone
tree, and the use of RAxML to estimate an ML tree on the
backbone alignment may help with the accuracy. If you obtain a
backbone tree in some other way, you should ensure that the
backbone tree and alignment are on the same exact set of taxa.
The backbone alignment is provided to \sepp using -a option, and should be in
the Fasta format. You can use refseq, available at
\url{http://www.ncbi.nlm.nih.gov/RefSeq/}, to convert between different alignment formats).
\item[Stats or info file:] this is the info file generated by RAxML (or phyml)
when it computes the ML tree on the backbone alignment (i.e., the
backbone tree).
This file is required by pplacer (run internally by \sepp)
in order to avoid re-estimation of ML parameters.
To be able to use \sepp you need to make sure you keep your
info file when you are generating the backbone tree.
If you do not have the info file (or if you used
some other software programs, such as PASTA~\cite{pasta-jcb}, to produce the
backbone tree), you can use RAxML's \ins{-f e} option to quickly
estimate the model parameters (including branch lengths) on your
backbone tree topology (see Section \ref{sec:backbone}).
The RAxML info file should be provided to \sepp using -r option.
\item[Fragments file:] this is a Fasta file containing the actual short fragments that are going to be placed. Fragments file should be given to \sepp using -f option.
\end{description}
\subsubsection{10\% rule}
Recall that \sepp operates by dividing the set of taxa
in the backbone alignment and tree into alignment subsets and placement subsets.
In the above run, we did not explicitly set the maximum subset sizes for alignment and placement subsets.
The choice of these two parameters affects accuracy on the one hand,
and computational resources on the other hand, as explained below.
\begin{itemize}
\item Decreasing alignment sizes should increase (and have increased in our experience) the accuracy of \sepp. On the other hand, smaller subsets increase the running time. This is because \sepp needs to score each fragment against all subsets independently, and therefore increasing the number of subsets adds to the running time. Note that extremely small subsets (i.e. less than 5 taxa) have not been tested extensively and are not recommended.
\item Increasing placement sizes should result in better
accuracy in general (although there could be exceptions).
If your placement tree is very large (thousands or
tens of thousands of leaves), the memory requirement of \pplacer,
and hence of \sepp, increases dramatically.
Reducing the placement size reduces the memory footprint, and hence
enables placement on larger trees given a fixed amount of memory available.
This would be one of the main motivations to reduce placement subset size.
Reducing the placement subset {\em can} result in reduced running time
as well, especially if your placement tree has thousands of taxa.
For smaller trees, the effect of the placement size on
the running time is not easily predicted, and is
practically of less interest.
\end{itemize}
By default, when alignment and placement subset sizes are
not explicitly specified by user, \sepp uses what we call the
``10\% rule" to automatically set those parameters. 10\% rule specifies that alignment and placement subset sizes should be both set to 10\% of the number of full length sequences (i.e. number of leaves in the backbone tree). The 10\% rule is just a heuristic setting we have found empirically to give a reasonable tradeoff {\em in general} between accuracy and computational requirements on the datasets we have tried. Users are encouraged to change subsets sizes based on their available computational resources, and the desired accuracy, according to the guidelines outlined above.
\subsubsection{Specifying subset sizes using -A and -P options}
Imagine we cannot wait 2 minutes to get results on our test dataset. We are going to increase the alignment subset so that \sepp runs faster. The test dataset included 65 full length sequences, and hence the 10\% rule amounts to alignment and placement subsets of maximum size 7. Since our toy example of 65 sequences is very small, it makes sense to increase the alignment size to a larger number (e.g. 10), and placement size to the entire dataset (i.e. 65). The maximum alignment and placement subset sizes are controlled with -A and -P options, respectively.
\begin{itemize}
\item Execute the following command to run \sepp with -A=10 and -P=65\\
\ins{run\_sepp.py -t mock/pyrg/sate.tre -r mock/pyrg/sate.tre.RAxML\_info -a mock/pyrg/sate.fasta -f mock/pyrg/pyrg.even.fas -A 10 -P 65 -o run2.A10P65}. \\
(\sepp should finish in less than a minute.)
\end{itemize}
In the above run note the -o option. This option controls the prefix of the output files generated by \sepp. We have not looked at \sepp output yet, and we will do so in a moment. For now, just be aware that \sepp generates a bunch of output files, and prefixes those with a given string, which defaults to ``output". These outputs are by default generated in the current directory, but that can be changed using the -d option. Had we not changed the output prefix, \sepp would have refused to run to avoid overwriting results of your previous run saved with ``output" prefix. Try this, and you would get the following error:\\
\ins{ Output directory [a\_drecotry] already contains files with prefix [output]...
Terminating to avoid loss of existing files.
}
%We are now ready to look the outputs of \sepp, but before we do that, let's start a \sepp run on a large dataset.
\subsection{Running \sepp on a larger dataset}
We are now going to run \sepp on a larger dataset.
Similar to the small ``pyrg'' dataset, our larger dataset is on fragments from a WGS sample of the HMMC mock community (\url{http://www.hmpdacc.org/HMMC/}),
but is based on a much larger marker gene called ``rpsS'' (obtained from \cite{metaphyler}).
The backbone alignment and tree are again estimated using \sate.
This dataset is available under \file{mock/rpsS} folder and has 1,277 full length sequences and 2101 fragments.
We are going to start a \sepp run on this dataset. Run the following command:\\
\ins{run\_sepp.py -t mock/rpsS/sate.tre -r mock/rpsS/sate.tre.RAxML\_info \\ -a mock/rpsS/sate.fasta -f mock/rpsS/rpsS.even.fas -o rpsS.out.default}\\
This run is going to take about 4 minutes. While this is running, we are going to look at \sepp outputs.
\section{\sepp output }
\sepp has two outputs: an alignment of fragments to full length sequences (or subsets of the full length sequences), and placements of fragments on the given backbone tree. When \sepp finishes running, it generates two or more output files, all prefixed with a string given using -o option, and placed in a directory given using -d option. One of these output files is \file{[prefix]\_placement.json} file. This is the results of the placement algorithm, and is in a format devised by the \pplacer software. This .json file is a human-readable text file. However in most cases you want to look at those results in a visualization tool. \pplacer comes with a suit of software tools called \guppy. \guppy reads a \file{.json} file, and among other things, produced nice visualizations of the results. We are going to first manually look at the plain \file{.json} file to understand its content, and then will use \guppy to visualize results.
\begin{itemize}
\item Open the file called \file{output\_placement.json} using a text editor.
\item Notice at the top of the file there is a newick tree with edges labeled with numbers inside brackets.
\item Following the tree, placement results
are given for each fragment. Everything between ``\{"and ``\}"
describes the placement of a single fragment. Each fragment can have multiple placements, with different likelihoods. Each line under the ``p" attribute indicates one placement of the fragment. The first value gives the edge label, the second value gives the log likelihood, the third value gives probability of that placement, and the final two values give the position on the edge where the fragment is placed, and the
length (as estimated using the maximum likelihood
calculation in pplacer)
of the pendant edge for the fragment.
\end{itemize}
Next we will use \guppy to turn this text file to a visualization of the results. Issue the following command to generate a tree that has fragments placed on the backbone tree based only on {\em the best placement} of each fragment.\\
\ins{\sbun/guppy tog --xml output\_placement.json}\\
This command generates a new file called \file{out\_placement.tog.xml}. This is an XML file in NexML format (\url{http://www.nexml.org/}) and can be opened with \arch software {available at \url{https://sites.google.com/site/cmzmasek/home/software/archaeopteryx}}.
%You can find \arch also on the USB drive provided to you by organizers of the workshop.
Run \arch and use \ins{File/Read Tree from File} to open \file{out\_placement.tog.xml}. Select ``colorize branch" checkbox on the right hand side panel (on newer versions of \arch, the option is ``visual styles/branch colors'' checkbox). The white branches represent the backbone tree, and branches in red correspond to the maximum likelihood placement of fragments on the backbone tree.
By now the run we started on the larger dataset (rpsS gene) should have finished as well (it takes about 5 minutes). Use \guppy to visualize the results of that run as well.
Please refer to \pplacer and \guppy documentation at \url{http://matsen.github.com/pplacer/generated_rst/pplacer.html} for more information about \file{.json} format and visualization options available in \guppy.
%Running SEPP using a configuration file
In addition to the placement file, an extended alignment file is also generated. This extended alignment shows how fragments are aligned to the backbone alignment. The extended alignment is a simple Fasta file, and can be viewed in any alignment visualization tool (e.g. JalView available at \url{http://www.jalview.org/}).
%Running SEPP on other datasets
\section{\sepp Obtaining Backbone Alignment and Trees}\label{sec:backbone}
In all examples given above, backbone alignment and trees were already estimated.
Our suggested way of obtaining backbone alignment and trees is through PASTA \cite{pasta-jcb},
which simultaneously estimates both an alignment and a tree based on unaligned full length sequences.
Please refer to \sate documentation for more information on running \sate.
In addition to an alignment and a tree (obtained from \sate or otherwise),
we also need to have a RAxML info file. If your backbone tree is estimated using
RAxML you already have the info file. Otherwise, you can optimize model parameters
on your backbone tree by running the following:\\\\
\ins{raxml -f e -t [backbone tree] -s [backbone alignment] -m GTRGAMMA -n [a name]}\\
This will optimize GTRGAMMA model parameters on your input alignment/tree pair and will
generate a info file (\file{RAxML\_info.a\_name}), that can be used with \sepp.
\section{\sepp Miscellaneous }
The following are some other points that are worth mentioning and testing.
\begin{itemize}
\item By default, the input is assumed to be DNA. To analyze amino acid datasets use \ins{-m} option (i.e. \ins{-m amino} for amino acid and \ins{-m rna} for RNA).
\item By default \sepp tries to use all available cores on your machine in each run. If you run multiple instances of \sepp simultaneously, or if you want it to use fewer cores, be sure
to set the number of cpus used by \sepp using \ins{-x} option. For example the following runs \sepp on the large dataset but with only 3 cpus used:\\
\ins{run\_sepp.py -t mock/rpsS/sate.tre -r mock/rpsS/sate.tre.RAxML\_info \\ -a mock/rpsS/sate.fasta -f mock/rpsS/rpsS.even.fas -x 3}\\
\item In addition to commandline, \sepp can be
controlled through a configuration file, passed to \sepp using \ins{-c} option. For example, to run \sepp using \file{config.run1} configuration file, use:\\
\ins{run\_sepp.py -c config.run1}\\
Commandline options can be specified in the configuration file under the \ins{[commandline]} section.
For specifying options under \ins{[commandline]}, you need to use their long format name (as show in the \sepp help invoked by \ins{-h}).
For example, to set input to ``rpsS'' dataset and the alignment size to 10 and number of cpus to 3 use:
\ins{
[commandline]\\
alignmentSize = 10 \\
tree= mock/rpsS/sate.tre \\
raxml = mock/rpsS/sate.tre.RAxML\_info \\
alignment = mock/rpsS/sate.fasta \\
fragment = mock/rpsS/rpsS.even.fas \\
cpu = 3}\\
Some extra information (not available in commandline) can also be configured in other sections of the configuration file. For example, \\
\ins{[pplacer]\\
path = /some/path\\
}
tells \sepp that \pplacer binaries can be found under \file{/some/path}.
An example config file is available as part of the distribution under the test directory (\file{test/unittest/data/simulated/sample.config}).
A main configuration file under \file{\{home\}/.sepp/main.config} is used
to store some basic configurations such as the location of extra programs, etc.
When conflicting options are given, precedence is with those provided through commandline, then those specified in config file provided using \ins{-c} option and finally those specified in the main config file. To test running from the config file, \ins{cd simulated} directory and run\\
\ins{run\_sepp.py -c test.config}
\item Currently \sepp has a built-in checkpointing functionality.
By default, this functionality is turned off, but can be turned on using \ins{-cp [checkpoint file name]},
and the time interval between two consecutive checkpoints can be adjusted using \ins{-cpi [checkpoint frequency in seconds]}.
Please note that checkpointing options have been tested only lightly and might have unknown issues.
\item For datasets containing both fragmentary and full-legnth sequences, the sequences much be separated out into a backbone set on the full-length sequences and a query set on the fragmentary sequences. From there, an alignment and tree can be estimated on the backbone set, and SEPP can be used to insert the query sequences back into the backbone tree. Run
\ins{split\_sequences.py -i mock/pyrg/data/mixed.fas -o split -t 100}\\
This will split the sequences into two files, split.full.fas (all sequences longer than 100 bps) and split.frag.fas *(all sequences shorter than or equal to 100 bps). From here, PASTA can be used to estimate the backbone alignment and RAxML can be used to estimate the backbone tree. The query sequences can be inserted into the tree.
\end{itemize}
\bibliographystyle{plain}
\bibliography{sepp}
\end{document}
|