phap - Phage Host Analysis Pipeline
A snakemake workflow that wraps various phage-host prediction tools.
When possible, tools and their dependencies are bundled in Singularity containers.
Current tools
Tool (source) | Publication/Preprint |
---|---|
RaFAh | Coutinho F. H. et al. 2020 |
vHuLK | Amgarten D. et al., 2020 |
VirHostMatcher-Net | Wang W. et al., 2020 |
WIsH | Galiez G. et al., 2017 |
Installation
Dependencies
To run the workflow your will need
-
snakemake > 5.x
(developed with5.30.1
) -
singularity >= 3.6
(developed with3.6.3
) -
biopython >= 1.78
(developed with1.78
)
Conda environemnt
It is recommended to use a conda environment.
The file environment.txt
can be used to recreate the complete environment
used during development.
The provided
environment.txt
contains an explicit list of all packages, produced withconda list -n hp --explicit > environment.txt
. This ensures all packages are exactly the same versions/builds, so we minimize the risk of running into dependencies issues
To get a working environment
# Clone this repo
$ git clone https://git.science.uu.nl/papanikos/phap.git
# Get in the dir
$ cd phap
# I am naming the environment `phap` here, you can call it whatever you like
# Note the long notation --file flag; -f will not work.
$ conda create -n phap --file=environment.txt
# Activate it - use the name you gave above, if it is different
$ conda activate phap
# The (hp) prefix shows we have activated it
# Check the snakemake version
(phap) $ snakemake --version
5.30.1
Configuration
Input data
The tools wrapped in this workflow expect phage sequences as input. You should try to make sure that the input sequences you want to analyze correspond to phage genomes/contigs (or at least viruses).
You can probably input any valid fasta file but the GIGO concept is probably applicable.
A separate workflow to identify phage/viral genomes/contigs is What the Phage.
The current workflow can handle multiple samples.
For each sample, all viral contigs to be analyzed should be provided as a
single multifasta (can be gz
ipped).
A mapping between sample ids and their corresponding fasta file is provided as
a samplesheet (see below).
Sample sheet
You must define a samplesheet with two comma (,
) separated columns and the
header sample,fasta
. Values from the sample
column must be unique and
are used as sample identifiers. Their corresponding fasta
values must be
valid paths to multifasta files with the phage sequences for that sample.
An example
$ cat samples.csv
sample,fasta
s01,/path/to/s01.fna
s02,/path/to/another.fna.gz
Note There is no need to follow any convention for the fasta file name to reflect the sample id. The values in the sample column are the ones to worry about, as these are the ones used as wildcards within the Snakefile.
You can
- Fill in the location of the samplesheet within the
config.yml
. - Drop the file in the workdir - Attention: It should be named
samples.csv
- Use
snakemake
's--config samplesheet=/path/to/my_samples.csv
when executing the wofkflow.
Models and data dependencies
-
RaFaH, vHULK For these tools there is no need to pre-download and setup anything - all data and software dependencies required for running them are bundled within the singularity image.
-
VirHostMatcher-Net, WIsH
Databases and models need to be downloaded from the VirHostMatcher data repo (see here). WIsH models for the 62,493 host genomes used in their paper are also provided and are used here for WIsH predictions.
Usage
Basic:
# From within this directory
# Make sure you have defined a samplesheet
(hp)$ snakemake --use-singularity -j16 \
--singularity-args "-B /path/to/databases/:/data"
where /path/to/database/
is the directory containing tables, WIsH models and
CRISPR blasts databases
Note
Binding the dir like this is required if the files are stored in some shared location and not on the local filesystem.
Output
All output is stored under a results
directory within the main workdir.
Results are stored per sample according to the sample ids you provided in the
sample sheet.
For each sample, results for each tool are stored in directories named after
the tool. An example looks like this:
results/A/
├── all_predictions.tsv
├── rafah
│ ├── A_CDS.faa
│ ├── A_CDS.fna
│ ├── A_CDS.gff
│ ├── A_CDSxMMSeqs_Clusters
│ ├── A_Genomes.fasta
│ ├── A_Genome_to_Domain_Score_Min_Score_50-Max_evalue_1e-05.tsv
│ ├── A_Ranger_Model_3_Predictions.tsv
│ ├── A_Seq_Info.tsv
│ └── predictions.tsv
├── tmp
│ ├── genomes
│ └── reflist.txt
├── vhmnet
│ ├── feature_values
│ ├── predictions
│ ├── predictions.tsv
│ └── tmp
├── vhulk
│ ├── predictions.tsv
│ └── results
└── wish
├── llikelihood.matrix
├── prediction.list
└── predictions.tsv
Per sample
-
all_predictions.tsv
: Contains the best prediction per contig (rows) for each tool along with its confidence/p-value/whatever single value each tool uses to evaluate its confidence in the prediction.
An example for three genomes:
contig vhulk_pred vhulk_score rafah_pred rafah_score vhmnet_pred vhmnet_score wish_pred wish_score
NC_005964.2 None 4.068828 Mycoplasma 0.461 Mycoplasma fermentans 0.9953 Bacteria;Tenericutes;Mollicutes;Mycoplasmatales;Mycoplasmataceae;Mycoplasma;Mycoplasma fermentans;Mycoplasma fermentans MF-I2 -1.2085700000000001
NC_015271.1 Escherichia_coli 1.0301523 Salmonella 0.495 Muricauda pacifica 0.9968 Bacteria;Proteobacteria;Gammaproteobacteria;Enterobacterales;Enterobacteriaceae;Raoultella;Raoultella sp. NCTC 9187;Raoultella sp. NCTC 9187 -1.3869200000000002
NC_023719.1 Bacillus 0.0012575098 Bacillus 0.55 Clostridium sp. LS 1.0000 Bacteria;Firmicutes;Clostridia;Clostridiales;Clostridiaceae;Clostridium;Clostridium beijerinckii;Clostridium beijerinckii -1.29454
-
tmp
directory- Contains one fasta file per input genome, along with other intermediate files necessary for a smooth execution of the workflow.
Per tool
-
rafah
- All files prefixed with
<sample_id>_
are the rafah's raw output -
predictions.tsv
: A selection of the 1st (Contig
) , 6th (Predicted_Host
) and 7th (Predicted_Host_Score
) columns from file<sample_id>_Seq_Info.tsv
- All files prefixed with
-
vhulk
-
results.csv
: Copy of theresults/sample/tmp/genomes/results/results.csv
-
predictions.tsv
: A selection of the 1st (BIN/genome
), 10th (final_prediction
) 11th (entropy
) columns from fileresults.csv
.
-
-
vhmnet
- Directories
feature_values
andpredictions
are the raw output - Directory
tmp
is a temporary dir written byVirHostMatcher-Net
for doing its magic. -
predictions.tsv
contain contig, host taxonomy and scores.
- Directories
-
wish
- Files
llikelihood.matrix
andprediction.list
are the raw output - File
predictions.tsv
has contig, host taxonomy and llikelihood scores.
- Files
Logs
Logs capturing stdout and stderr during execution of each rule can be found in
workdir/logs/<sample_id>/*.log
files.