Variety of transformations for abundance data, stored in assay
.
See details for options.
transformSamples(
x,
assay.type = "counts",
assay_name = NULL,
method = c("alr", "chi.square", "clr", "frequency", "hellinger", "log", "log10",
"log2", "normalize", "pa", "rank", "rclr", "relabundance", "rrank", "total"),
name = method,
...
)
# S4 method for SummarizedExperiment
transformSamples(
x,
assay.type = "counts",
assay_name = NULL,
method = c("alr", "chi.square", "clr", "frequency", "hellinger", "log", "log10",
"log2", "normalize", "pa", "rank", "rclr", "relabundance", "rrank", "total"),
name = method,
pseudocount = FALSE,
...
)
transformAssay(
x,
assay.type = "counts",
assay_name = NULL,
method = c("alr", "chi.square", "clr", "frequency", "hellinger", "log", "log10",
"log2", "max", "normalize", "pa", "range", "rank", "rclr", "relabundance", "rrank",
"standardize", "total", "z"),
MARGIN = "samples",
name = method,
pseudocount = FALSE,
...
)
transformCounts(x, ...)
# S4 method for SummarizedExperiment
transformAssay(
x,
assay.type = "counts",
assay_name = NULL,
method = c("alr", "chi.square", "clr", "frequency", "hellinger", "log", "log10",
"log2", "max", "normalize", "pa", "range", "rank", "rclr", "relabundance", "rrank",
"standardize", "total", "z"),
MARGIN = "samples",
name = method,
pseudocount = FALSE,
...
)
transformFeatures(
x,
assay.type = "counts",
assay_name = NULL,
method = c("frequency", "log", "log10", "log2", "max", "pa", "range", "standardize",
"z"),
name = method,
pseudocount = FALSE,
...
)
# S4 method for SummarizedExperiment
transformFeatures(
x,
assay.type = "counts",
assay_name = NULL,
method = c("frequency", "log", "log10", "log2", "max", "pa", "range", "standardize",
"z"),
name = method,
pseudocount = FALSE,
...
)
ZTransform(x, MARGIN = "features", ...)
# S4 method for SummarizedExperiment
ZTransform(x, MARGIN = "features", ...)
relAbundanceCounts(x, ...)
# S4 method for SummarizedExperiment
relAbundanceCounts(x, ...)
A
SummarizedExperiment
object.
A single character value for selecting the
assay
to be
transformed.
a single character
value for specifying which
assay to use for calculation.
(Please use assay.type
instead. At some point assay_name
will be disabled.)
A single character value for selecting the transformation method.
A single character value specifying the name of transformed abundance table.
additional arguments passed on to vegan:decostand
:
ref_vals
: A single value which will be used to fill
reference sample's column in returned assay when calculating alr.
(default: ref_vals = NA
)
TRUE, FALSE, or a numeric value. When TRUE,
automatically adds the minimum positive value of assay.type
.
When FALSE, does not add any pseudocount (pseudocount = 0).
Alternatively, a user-specified numeric value can be added as pseudocount.
A single character value for specifying whether the
transformation is applied sample (column) or feature (row) wise.
(By default: MARGIN = "samples"
)
transformAssay
returns the input object x
, with a new
transformed abundance table named name
added in the assay
.
These transformCount
function provides a variety of options for transforming abundance data.
The transformed data is calculated and stored in a new assay
. The previously available
wrappers transformSamples, transformFeatures
ZTransform, and relAbundanceCounts have been deprecated.
The transformAssay
provides sample-wise (column-wise) or feature-wise
(row-wise) transformation to the abundance table
(assay) based on specified MARGIN
.
The available transformation methods include:
'alr' Additive log ratio (alr) transformation, please refer to
decostand
for details.
'chi.square' Chi square transformation, please refer to
decostand
for details.
'clr' Centered log ratio (clr) transformation, please refer to
decostand
for details.
'frequency' Frequency transformation, please refer to
decostand
for details.
'hellinger' Hellinger transformation, please refer to
decostand
for details.
'log' Logarithmic transformation, please refer to
decostand
for details.
'log10' log10 transformation can be used for reducing the skewness of the data. $$log10 = \log_{10} x$$ where \(x\) is a single value of data.
'log2' log2 transformation can be used for reducing the skewness of the data. $$log2 = \log_{2} x$$ where \(x\) is a single value of data.
'normalize' Make margin sum of squares equal to one. Please refer to
decostand
for details.
'pa' Transforms table to presence/absence table. Please refer to
decostand
for details.
'rank' Rank transformation, please refer to
decostand
for details.
'rclr' Robust clr transformation, please refer to
decostand
for details.
'relabundance' Relative transformation (alias for 'total'), please refer to
decostand
for details.
'rrank' Relative rank transformation, please refer to
decostand
for details.
'standardize' Scale 'x' to zero mean and unit variance (alias for
'z'), please refer to decostand
for details.
'total' Divide by margin total (alias for
'relabundance'), please refer to
decostand
for details.
'z' Z transformation (alias for 'standardize'),
please refer to decostand
for details.
data(esophagus, package="mia")
tse <- esophagus
# By specifying 'method', it is possible to apply different transformations,
# e.g. compositional transformation.
tse <- transformAssay(tse, method = "relabundance")
# The target of transformation can be specified with "assay.type"
# Pseudocount can be added by specifying 'pseudocount'.
# Perform CLR with smallest positive value as pseudocount
tse <- transformAssay(tse, assay.type = "relabundance", method = "clr",
pseudocount = TRUE
)
#> A pseudocount of 0.00392156862745098 was applied.
head(assay(tse, "clr"))
#> B C D
#> 59_8_22 3.4161459 2.0677373 1.10139252
#> 59_5_13 -0.7397308 0.1706173 -0.81874812
#> 59_8_12 0.5165446 1.7110624 -0.81874812
#> 65_3_22 -0.7397308 0.1706173 -0.04661252
#> 65_5_1 -0.7397308 -0.9279950 -0.04661252
#> 65_1_10 -0.7397308 -0.9279950 -0.04661252
# With MARGIN, you can specify the if transformation is done for samples or
# for features. Here Z-transformation is done feature-wise.
tse <- transformAssay(tse, method = "z", MARGIN = "features")
#> Warning: result contains NaN, perhaps due to impossible mathematical
#> operation
head(assay(tse, "z"))
#> B C D
#> 59_8_22 1.1000638 -0.2460669 -0.8539969
#> 59_5_13 -0.5773503 1.1547005 -0.5773503
#> 59_8_12 -0.4285714 1.1428571 -0.7142857
#> 65_3_22 -1.0000000 1.0000000 0.0000000
#> 65_5_1 -0.5773503 -0.5773503 1.1547005
#> 65_1_10 -0.5773503 -0.5773503 1.1547005
# Name of the stored table can be specified.
tse <- transformAssay(tse, method="hellinger", name="test")
head(assay(tse, "test"))
#> B C D
#> 59_8_22 0.49629167 0.27296484 0.15109947
#> 59_5_13 0.00000000 0.08856149 0.00000000
#> 59_8_12 0.09925833 0.22578838 0.00000000
#> 65_3_22 0.00000000 0.08856149 0.06757374
#> 65_5_1 0.00000000 0.00000000 0.06757374
#> 65_1_10 0.00000000 0.00000000 0.06757374
# pa returns presence absence table.
tse <- transformAssay(tse, method = "pa")
head(assay(tse, "pa"))
#> B C D
#> 59_8_22 1 1 1
#> 59_5_13 0 1 0
#> 59_8_12 1 1 0
#> 65_3_22 0 1 1
#> 65_5_1 0 0 1
#> 65_1_10 0 0 1
# rank returns ranks of taxa.
tse <- transformAssay(tse, method = "rank")
head(assay(tse, "rank"))
#> B C D
#> 59_8_22 27.0 30 33.0
#> 59_5_13 0.0 13 0.0
#> 59_8_12 16.5 28 0.0
#> 65_3_22 0.0 13 11.5
#> 65_5_1 0.0 0 11.5
#> 65_1_10 0.0 0 11.5
# In order to use other ranking variants, modify the chosen assay directly:
assay(tse, "rank_average", withDimnames = FALSE) <- colRanks(assay(tse, "counts"),
ties.method="average",
preserveShape = TRUE)