Detect differential ChIP signal using one of many approaches
fitAssayDiff(x, ...)
# S4 method for class 'SummarizedExperiment'
fitAssayDiff(
x,
assay = "counts",
design = NULL,
coef = NULL,
lib.size = "totals",
method = c("qlf", "lt", "wald"),
norm = c("none", "TMM", "RLE", "TMMwsp", "upperquartile"),
groups = NULL,
fc = 1,
lfc = log2(fc),
asRanges = FALSE,
offset = NULL,
weighted = FALSE,
...,
null = c("interval", "worst.case"),
robust = FALSE,
type = c("apeglm", "ashr", "normal")
)a SummarizedExperiment object
Passed to normLibSizes and
The assay to use for analysis
The design matrix to use for analysis
The required column from the design matrix
The column within the colData element which contains the library size information. If set to NULL, column summaries will be used.
the analytic method to be used. Can be 'qlf' which will fit counts using the glmQLFit strategy , or 'lt' which fits the limma-trend model on logCPM, or pre-processed logCPM values. Setting method = 'wald' will call nbinomWaldTest
The normalisation strategy to use when running the glmQLF model or the Wald test. The value 'none' relies solely on library-size normalisation, and is the default. All methods available in normLibSizes are implemented. Ignored when using method = "lt"
character(1) If a column name is supplied here, group-based normalisation will be applied to GLM models treating data in this column as a grouping factor. Ignored when using method = "lt"
logical(1). By default, the returned object will be a
SummarizedExperiment object with the results added to the rowData
element. Setting asRanges = TRUE will only return the GRanges object from
this element
If provided will be used as the offset when a DGEList object is created during model fitting for method = 'qlf'
logical(1) Passed to normLibSizes. Only used when applying a TMM-type normalisation strategy
Passed to glmTreat glmQLFit when method = "qlf". If method = "lt", instead passed to lmFit
Passed to lfcShrink
A SummarizedExperiment object with results set as the rowData element.
Any existing columns not contained in the differential ChIP results will be
retained.
Results from testing will contain logCPM, logFC, PValue and any t/F
statistic as appropriate, along with an FDR-adjusted p-value.
If specifying a range-based H0 by setting lfc != 0, an additional column p_mu0 will be included which is the p-value for the point H0: logFC = 0. These are not used for FDR-adjusted p-values but can be helpful when integrating multiple ChIP targets due to the increase in false-negatives when using a range-based H0, and when requiring more accurate identification of truly unchanged sites, as opposed to those which simply fail to achieve significance using a range-based H0 where arbitrary cutoff values are used.
Starting with a SummarizedExperiment object this function fits either a glmQLFit model to count data, performs the nbinomWaldTest on count data, or applies the limma-trend model to logCPM data.
If fitting Generalised Linear Models via glmQLFit, options for normalisation are "none", which normalises to library size. Existing library sizes are commonly found in the "totals" column of the colData element and this is attempted by default. All methods provided in normLibSizes are also implemented, with the added possibility of normalising within groups instead of across the entire dataset. To enable this, the column with the grouping factor is expected to be in the colData element and is simply called by column name. No normalisation is applied when using the limma-trend model, as this allows for previous normalisation strategies to be performed on the data.
If testing with nbinomWaldTest, applying RLE normalisation without groups, and using colSums for library sizes (instead of total alignments), the standard normalisation factors from estimateSizeFactors will be used. In all other scenarios, normalisation factors as returned by normLibSizes will be used. The fitType is set to 'local' when estimating dispersions, and this can be easily modified by passing fitType via the dot arguments. Results are additionally returned after applying lfcShrink, including the svalue returned by this approach.
Normalising to ChIP Input samples is not yet implemented. Similarly, the use of offsets when applying the Wald test is not yet implemented.
Range-based hypothesis testing is implemented using glmTreat or treat. Setting fc to 1 (or lfc to 0) will default to a point-based null hypothesis, equivalent to either glmQLFTest (method = "qlf") or eBayes (method = "lt"). When applying nbinomWaldTest, lfcShrink will be applied.
It should also be noted that this is primarily a convenience function and if requiring intermediate output from any steps, then these can be run individually as conventionally specified.
nrows <- 200; ncols <- 6
counts <- matrix(runif(nrows * ncols, 1, 1e4), nrows)
colnames(counts) <- paste0("Sample_", seq_len(ncols))
df <- DataFrame(treat = c("A", "A", "A", "B", "B", "B"))
df$treat <- as.factor(df$treat)
se <- SummarizedExperiment(
assays = SimpleList(counts = counts), colData = df
)
X <- model.matrix(~treat, colData(se))
se <- fitAssayDiff(se, design = X, lib.size = NULL)
#> Creating DGE list...
#> Calculating experiment-wide normalisation factors...
#> Estimating dispersions...
#> Running glmQLFit...
rowData(se)
#> DataFrame with 200 rows and 5 columns
#> logFC logCPM F PValue FDR
#> <numeric> <numeric> <numeric> <numeric> <numeric>
#> 1 4.3118042 11.7021 10.79729617 0.00272349 0.544698
#> 2 0.0832925 12.2407 0.00682725 0.93473200 0.996024
#> 3 1.6788140 11.6748 2.04323371 0.16387523 0.996024
#> 4 -1.1910443 12.0524 1.32867823 0.25872306 0.996024
#> 5 0.5183814 12.0442 0.24296590 0.62589044 0.996024
#> ... ... ... ... ... ...
#> 196 0.0906305 12.7687 0.02760966 0.869216 0.996024
#> 197 1.0957446 12.3006 1.21717268 0.279247 0.996024
#> 198 -1.1430050 12.6072 2.87907802 0.100755 0.996024
#> 199 0.1462205 12.3253 0.02519284 0.875019 0.996024
#> 200 0.0459551 11.8591 0.00124162 0.972139 0.996024