Skip to content

Hidden state prediction

Gavin Douglas edited this page May 22, 2018 · 22 revisions

PICRUSt2 wraps the castor R package to run hidden-state prediction (hsp) to predict gene family abundances. See below how the --chunk_size and -p options can be used to maximize memory and run-time efficiency.

Hidden-state prediction for 16S copy number and E.C. numbers can be run with these commands:

hsp.py -i 16S -t OUTPUT.tre -o 16S_predicted -p 10 -n

hsp.py -i EC -t OUTPUT.tre -o EC_predicted -p 10

An important step is to look at the NSTI values included in the 16S output table (which was indicated by the -n option). Sequences with extremely high NSTI values (e.g. > 1) should be removed since they will only add noise.

The input arguments and options are:

  • -t TREEFILE: Newick tree with study sequences placed amongst reference sequences.
  • -i TRAIT_OPTION: Which default pre-calculated count table to use (one of '16S', 'COG', 'EC', 'KO', 'PFAM', 'TIGRFAM')
  • -o PREFIX: Prefix for output files: RDS (R object containing state probabilities), predicted counts, and optionally a table of CIs.
  • --observed_trait_table TRAIT_COUNTFILE: Trait file to use if a non-default file is needed (most users should use one of the default trait options above).
  • -m METHOD: Hidden-state prediction method to use, which needs to be one of: maximum parsimony (mp), empirical probabilities (emp_prob), Markov model (mk_model), subtree averaging (subtree_average), phylogenetic independent contrast (pic), or squared-change parsimony (scq).
  • -p INT: Number of processes to run in parallel.
  • --chunk_size INT: Number of gene families to read in for a given processor. Note that increasing this value will not always speed up execution! The trait table is split into subsets of size chunk_size before running hidden state prediction to reduce memory usage. Each subset can be run through on a different processor, but the number of simultaneous processors is equal to the number of data subsets! For example, if there are 1500 total gene families and the chunk_size is 500, the maximum number of processes that can be run is parallel is 3 (which will be the case even if you set something higher!).
  • -n: Indicates that Nearest-sequenced taxon index (NSTI) values should be calculated. This metric can be used to identify study sequences that are highly distant from all reference sequences (the predictions for these sequences are less reliable!). There is no clear cut-off for a high NSTI values, but a good rule of thumb is that sequences placed with NSTI > 0.15 will be less reliable.

Clone this wiki locally