testNhoods | R Documentation |
This will perform differential neighbourhood abundance testing after cell counting.
x |
A |
design |
A |
design.df |
A |
min.mean |
A scalar used to threshold neighbourhoods on the minimum average cell counts across samples. |
model.contrasts |
A string vector that defines the contrasts used to perform DA testing. |
fdr.weighting |
The spatial FDR weighting scheme to use. Choice from max,
neighbour-distance, graph-overlap or k-distance (default). If |
robust |
If robust=TRUE then this is passed to edgeR and limma which use a robust
estimation for the global quasilikelihood dispersion distribution. See |
norm.method |
A character scalar, either |
reduced.dim |
A character scalar referring to the reduced dimensional slot used to compute distances for the spatial FDR. This should be the same as used for graph building. |
This function wraps up several steps of differential abundance testing using
the edgeR
functions. These could be performed separately for users
who want to exercise more contol over their DA testing. By default this
function sets the lib.sizes
to the colSums(x), and uses the
Quasi-Likelihood F-test in glmQLFTest
for DA testing. FDR correction
is performed separately as the default multiple-testing correction is
inappropriate for neighbourhoods with overlapping cells.
A data.frame
of model results, which contain:
logFC
:Numeric, the log fold change between conditions, or for an ordered/continous variable the per-unit change in (normalized) cell counts per unit-change in experimental variable.
logCPM
:Numeric, the log counts per million (CPM), which equates to the average log normalized cell counts across all samples.
F
:Numeric, the F-test statistic from the quali-likelihood F-test
implemented in edgeR
.
PValue
:Numeric, the unadjusted p-value from the quasi-likelihood F-test.
FDR
:Numeric, the Benjamini & Hochberg false discovery weight
computed from p.adjust
.
Nhood
:Numeric, a unique identifier corresponding to the specific graph neighbourhood.
SpatialFDR
:Numeric, the weighted FDR, computed to adjust for spatial graph overlaps between neighbourhoods. For details see graphSpatialFDR.
Mike Morgan
library(SingleCellExperiment) ux.1 <- matrix(rpois(12000, 5), ncol=400) ux.2 <- matrix(rpois(12000, 4), ncol=400) ux <- rbind(ux.1, ux.2) vx <- log2(ux + 1) pca <- prcomp(t(vx)) sce <- SingleCellExperiment(assays=list(counts=ux, logcounts=vx), reducedDims=SimpleList(PCA=pca$x)) milo <- Milo(sce) milo <- buildGraph(milo, k=20, d=10, transposed=TRUE) milo <- makeNhoods(milo, k=20, d=10, prop=0.3) milo <- calcNhoodDistance(milo, d=10) cond <- rep("A", ncol(milo)) cond.a <- sample(1:ncol(milo), size=floor(ncol(milo)*0.25)) cond.b <- setdiff(1:ncol(milo), cond.a) cond[cond.b] <- "B" meta.df <- data.frame(Condition=cond, Replicate=c(rep("R1", 132), rep("R2", 132), rep("R3", 136))) meta.df$SampID <- paste(meta.df$Condition, meta.df$Replicate, sep="_") milo <- countCells(milo, meta.data=meta.df, samples="SampID") test.meta <- data.frame("Condition"=c(rep("A", 3), rep("B", 3)), "Replicate"=rep(c("R1", "R2", "R3"), 2)) test.meta$Sample <- paste(test.meta$Condition, test.meta$Replicate, sep="_") rownames(test.meta) <- test.meta$Sample da.res <- testNhoods(milo, design=~Condition, design.df=test.meta[colnames(nhoodCounts(milo)), ], norm.method="TMM") da.res
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.