philr {philr} | R Documentation |
This is the main function for building the phylogenetic ILR basis, calculating the weightings (of the parts and the ILR coordinates) and then transforming the data.
philr( x, tree = NULL, sbp = NULL, part.weights = "uniform", ilr.weights = "uniform", return.all = FALSE, abund_values = "counts" )
x |
matrix of data to be transformed (samples are rows, compositional parts are columns) - zero must be dealt with either with pseudocount, multiplicative replacement, or another method. |
tree |
a |
sbp |
(Optional) give a precomputed sbp matrix |
part.weights |
weightings for parts, can be a named vector with names
corresponding to
|
ilr.weights |
weightings for the ILR coordiantes can be a named vector
with names corresponding to names of internal nodes of
|
return.all |
return all computed parts (e.g., computed sign
matrix( |
abund_values |
A single character value for selecting the
|
This is a utility function that pulls together a number of other functions
in philr
. The steps that are executed are as follows:
Create sbp (sign matrix) if not given
Create parts weightings if not given
Shift the dataset with respect to the new reference measure (e.g., part weightings)
Create the basis contrast matrix from the sign matrix and the reference measure
Transform the data based on the contrast matrix and the reference measure
Calculate the specified ILR weightings and multiply each balance by the corresponding weighting
Note for both the reference measure (part weightings) and the ILR weightings,
specifying 'uniform'
will give the same results as not weighting at
all.
Note that some of the prespecified part.weights assume x
is given as
counts and not as relative abundances. Except in this case counts or relative
abundances can be given.
The tree argument is ignored if the x
argument is
assay
or
assay
These objects can include a phylogenetic tree. If the phylogenetic tree
is missing from these objects, it should be integrated directly in these
data objects before running philr
. Alternatively, you can always
provide the abundance matrix and tree separately in their standard formats.
If you have a
assay
,
this can be converted into
assay
,
to incorporate tree information.
matrix if return.all=FALSE
, if return.all=TRUE
then
a list is returned (see above).
Justin Silverman; S3 methods by Leo Lahti
# Prepare example data tr <- named_rtree(5) x <- t(rmultinom(10,100,c(.1,.6,.2,.3,.2))) + 0.65 # add a small pseudocount colnames(x) <- tr$tip.label philr(x, tr, part.weights='enorm.x.gm.counts', ilr.weights='blw.sqrt', return.all=FALSE) # Running philr on a TreeSummarizedExperiment object ## Prepare example data library(mia) library(tidyr) data(GlobalPatterns, package="mia") ## Select prevalent taxa tse <- GlobalPatterns %>% subsetByPrevalentTaxa( detection = 3, prevalence = 20/100, as_relative = FALSE) ## Pick taxa that have notable abundance variation across sammples variable.taxa <- apply(assay(tse, "counts"), 1, function(x) sd(x)/mean(x) > 3.0) tse <- tse[variable.taxa,] # Collapse the tree tree <- ape::keep.tip(phy = rowTree(tse), tip = rowLinks(tse)$nodeNum) rowTree(tse) <- tree ## Add a new assay with a pseudocount assays(tse)$counts.shifted <- assay(tse, "counts") + 1 ## Run philr for TreeSummarizedExperiment object ## using the pseudocount data res.tse <- philr(tse, part.weights='enorm.x.gm.counts', ilr.weights='blw.sqrt', return.all=FALSE, abund_values="counts.shifted") # Running philr on a phyloseq object pseq <- makePhyloseqFromTreeSummarizedExperiment(tse) res.pseq <- philr(pseq, part.weights='enorm.x.gm.counts', ilr.weights='blw.sqrt', return.all=FALSE)