ScpModel-DifferentialAnalysis | R Documentation |
Differential abundance analysis assess the statistical significance of the differences observed between group of samples of interest.
scpDifferentialAnalysis(object, coefficients = NULL, contrasts = NULL, name)
scpDifferentialAggregate(differentialList, fcol, ...)
scpVolcanoPlot(
differentialList,
fdrLine = 0.05,
top = 10,
by = "padj",
decreasing = FALSE,
textBy = "feature",
pointParams = list(),
labelParams = list()
)
object |
An object that inherits from the
|
coefficients |
A |
contrasts |
A |
name |
A |
differentialList |
A list of tables returned by
|
fcol |
A |
... |
Further arguments passed to
|
fdrLine |
A |
top |
A |
by |
A |
decreasing |
A |
textBy |
A |
pointParams |
A |
labelParams |
A |
scpDifferentialAnalysis()
performs a t-test on a coefficient of
interest or on coefficients that distinguish two groups of cells
of interests (provided as a contrast). Contrasts must be provided
as a list where each element is a three-element character vector.
The first element of the vector provides the name of the
categorical variable to test, the second element provides the
name of the first category (that is one of the factor levels), the
third element provides the name of the other category to compare.
Numerical variables can be tested by providing the coefficient
argument, that is the name of the coefficient associated to that
numerical variable. The statistical tests are conducted for each
feature independently. The p-values are adjusted using
IHW::ihw()
, where each test is weighted using the feature
intercept (that is the average baseline intensity). The function
returns a list of DataFrame
s with one table for each test
contrast and/or coefficient. It provides the adjusted p-values and
the estimates. For contrast, the estimates represent the estimated
log fold changes between the groups. For coefficients, the
estimates are the estimated slopes.
scpDifferentialAggregate()
combines the differential abundance
analysis results for groups of features. This is useful, for
example, to return protein-level results when data is modelled at
the peptide level. The function heavily relies on the approaches
implemented in metapod::combineGroupedPValues()
. The p-values
are combined into a single value using one of the following
methods: Simes' method
(default), Fisher's method, Berger's method, Pearson's method,
minimum Holm's approach, Stouffer's Z-score method, and
Wilkinson's method. We refer to the metapod
documentation for
more details on the assumptions underlying each approach. The
estimates are combined using the representative estimate, as
defined by metapod
. Which estimate is representative depends on
the selected combination method. The function takes the list of
tables generated by scpDifferentialAnalysis()
and returns a new
list of DataFrame
s with aggregated results.
scpAnnotateResults()
adds annotations to the differential abundance
analysis results. The annotations are added to all elements of the
list returned by ()
. See the associated
man page for more information.
scpVolcanoPlot()
takes the list of tables generated by
scpDifferentialAnalysis()
and returns a ggplot2
scatter plot.
The plots show the adjusted p-values with respect to the estimate.
A horizontal bar also highlights the significance threshold
(defaults to 5%, fdrLine
). The top (default 10) features with lowest
p-values are labeled on the plot. You can control which features
are labelled using the top
, by
and decreasing
arguments.
Finally, you can change the point and label aesthetics thanks to
the pointParams
and the labelParams
arguments, respectively.
Christophe Vanderaa, Laurent Gatto
ScpModel-Workflow to run a model on SCP data upstream of differential abundance analysis.
scpAnnotateResults()
to annotate analysis of variance results.
library("patchwork")
library("ggplot2")
data("leduc_minimal")
## Add n/p ratio information in rowData
rowData(leduc_minimal)$npRatio <-
scpModelFilterNPRatio(leduc_minimal, filtered = FALSE)
####---- Run differential abundance analysis ----####
(res <- scpDifferentialAnalysis(
leduc_minimal, coefficients = "MedianIntensity",
contrasts = list(c("SampleType", "Melanoma", "Monocyte"))
))
## IHW return a message because of the example data set has only few
## peptides, real dataset should not have that problem.
####---- Annotate results ----####
## Add peptide annotations available from the rowData
res <- scpAnnotateResults(
res, rowData(leduc_minimal),
by = "feature", by2 = "Sequence"
)
####---- Plot results ----####
scpVolcanoPlot(res, textBy = "gene") |>
wrap_plots(guides = "collect")
## Modify point and label aesthetics
scpVolcanoPlot(
res, textBy = "gene", top = 20,
pointParams = list(aes(colour = npRatio), alpha = 0.5),
labelParams = list(size = 2, max.overlaps = 20)) |>
wrap_plots(guides = "collect")
####---- Aggregate results ----####
## Aggregate to protein-level results
byProteinDA <- scpDifferentialAggregate(
res, fcol = "Leading.razor.protein.id"
)
scpVolcanoPlot(byProteinDA) |>
wrap_plots(guides = "collect")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.