plotAbundanceDensity.Rd
This function plots abundance of the most abundant taxa.
plotAbundanceDensity(object, ...)
# S4 method for SummarizedExperiment
plotAbundanceDensity(
object,
layout = c("jitter", "density", "point"),
assay.type = assay_name,
assay_name = "counts",
n = min(nrow(object), 25L),
colour_by = NULL,
shape_by = NULL,
size_by = NULL,
order_descending = TRUE,
...
)
a
SummarizedExperiment
object.
additional parameters for plotting.
xlab a single character value for selecting the x-axis label.
(default: xlab = assay.type
)
ylab a single character value for selecting the y-axis label.
ylab
is disabled when layout = "density"
.
(default: ylab = "Taxa")
point_alpha a numeric value from range 0 to 1 selecting the transparency of
colour in jitter
and point
plot. (default: point_alpha = 0.6
)
point_shape a positive integer value selecting the shape of point in
jitter
and point
plot. (default: point_shape = 21
)
point_size a positive numeric value selecting the size of point in
jitter
and point
plot. (default: point_size = 2
)
add_legend a boolean value selecting if legend is added.
(default: add_legend = TRUE
)
flipped a boolean value selecting if the orientation of plot is changed
so that x-axis and y-axis are swapped. (default flipped = FALSE
)
add_x_text a boolean value selecting if text that represents values is included
in x-axis. (default: add_x_text = TRUE
)
See mia-plot-args
for more details i.e. call help("mia-plot-args")
a single character value for selecting the layout of the plot.
There are three different options: jitter
, density
, and
point
plot. (default: layout = "jitter"
)
a single character value for selecting the
assay
to be
plotted. (default: assay.type = "counts"
)
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 positive integer specifying the number of the most abundant taxa
to show. (default: n = min(nrow(object), 25L)
)
a single character value defining a column from
colData
, that is used to color plot. Must be a value of
colData()
function. (default: colour_by = NULL
)
a single character value defining a column from
colData
, that is used to group observations to different point shape
groups. Must be a value of colData()
function. shape_by
is
disabled when layout = "density"
. (default: shape_by = NULL
)
a single character value defining a column from
colData
, that is used to group observations to different point size
groups. Must be a value of colData()
function. size_by
is
disabled when layout = "density"
. (default: size_by = NULL
)
TRUE
, FALSE
orNA
: Should the
results be ordered in a descending order? If NA
is given the order
as found in object
for the n
most abundant taxa is used.
(default: order_descending = TRUE
)
A ggplot2
object
This function plots abundance of the most abundant taxa. Abundance can be plotted as a jitter plot, a density plot, or a point plot. By default, x-axis represents abundance and y-axis taxa. In a jitter and point plot, each point represents abundance of individual taxa in individual sample. Most common abundances are shown as a higher density.
A density plot can be seen as a smoothened bar plot. It visualized distribution of abundances where peaks represent most common abundances.
tse <- microbiomeDataSets::atlas1006()
#> see ?microbiomeDataSets and browseVignettes('microbiomeDataSets') for documentation
#> loading from cache
#> see ?microbiomeDataSets and browseVignettes('microbiomeDataSets') for documentation
#> loading from cache
#> see ?microbiomeDataSets and browseVignettes('microbiomeDataSets') for documentation
#> loading from cache
# Plots the abundances of 25 most abundant taxa. Jitter plot is the default option.
plotAbundanceDensity(tse, assay.type = "counts")
# Counts relative abundances
tse <- transformAssay(tse, method = "relabundance")
# Plots the relative abundance of 10 most abundant taxa.
# "nationality" information is used to color the points. X-axis is log-scaled.
plotAbundanceDensity(tse, layout = "jitter", assay.type = "relabundance",
n = 10, colour_by = "nationality") +
scale_x_log10()
# Plots the relative abundance of 10 most abundant taxa as a density plot.
# X-axis is log-scaled
plotAbundanceDensity(tse, layout = "density", assay.type = "relabundance",
n = 10 ) +
scale_x_log10()
# Plots the relative abundance of 10 most abundant taxa as a point plot.
# Point shape is changed from default (21) to 41.
plotAbundanceDensity(tse, layout = "point", assay.type = "relabundance", n = 10,
point_shape = 41)
# Plots the relative abundance of 10 most abundant taxa as a point plot.
# In addition to colour, groups can be visualized by size and shape in point plots,
# and adjusted for point size
plotAbundanceDensity(tse, layout = "point", assay.type = "relabundance", n = 10,
shape_by = "sex", size_by = "time", point_size=1)
#> Warning: Removed 370 rows containing missing values or values outside the scale range
#> (`geom_point()`).
# Ordering via order_descending
plotAbundanceDensity(tse, assay.type = "relabundance",
order_descending = FALSE)
# for custom ordering set order_descending = NA and order the input object
# to your wishes
plotAbundanceDensity(tse, assay.type = "relabundance",
order_descending = NA)