Create design & contrast matrix for differential state analysis

[Stein-ErikG]

I have a mass cytometry dataset of samples collected from patients with AML (n=30). There is one sample from each patient, and the patients can be divided into two groups based on initial response to experimental therapy (21 responders vs 9 non-responders). I have been using the CATALYST pipeline and the diffcyt() function to do differential state analysis. Using this I find lot if very significant differences between the two groups. But when I use a different script that just do a lot of two sided T-tests (and then FDR adjustment), I find nothing significant.

When I use diffcyt() and my t-test script on a different and paired dataset, they find similar significant differences. So, I surmise that my design and/or contrast is not correctly set up? Can anyone help me? I hope my R output and code snippet below is informative! 😊

design <- createDesign(ei(sce), cols_design = “condition”)
(Intercept) conditionRD
1 1 1
2 1 1
3 1 0
4 1 1
5 1 0
6 1 0
7 1 0
8 1 0
9 1 0
10 1 0
11 1 0
12 1 0
13 1 0
14 1 0
15 1 0
16 1 1
17 1 1
18 1 1
19 1 0
20 1 1
21 1 0
22 1 1
23 1 1
24 1 0
25 1 0
26 1 0
27 1 0
28 1 0
29 1 0
30 1 0
attr(,"assign")
[1] 0 1
attr(,"contrasts")
attr(,"contrasts")$condition
[1] "contr.treatment"

contrast <- createContrast(c(0,1))
[,1]
[1,] 0
[2,] 1

I am using the diffcyt() (https://rdrr.io/bioc/diffcyt/man/diffcyt.html) function from the diffcyt pipeline, with mostly default conditions. Ive pasted the code im using below:

design <- createDesignMatrix(ei(sce), cols_design = cond)
contrast <- createContrast(c(0, 1))

res_DS <- diffcyt(sce,
clustering_to_use = clust,
analysis_type = "DS",
method_DS = "diffcyt-DS-limma",
design = design,
contrast = contrast,
verbose = FALSE)
tbl_DS <- rowData(res_DS$res)

plotDiffHeatmap(sce, tbl_DS, all = TRUE, top_n = 30, fdr = 0.1,
sort_by = "padj", col_anno = cond, normalize =TRUE)

I assume that “diffcyt-DS-limma” uses limma::lmFit. It looks to me that diffcyt() (https://rdrr.io/bioc/diffcyt/src/R/diffcyt_wrapper.R) uses testDS_limma.R with my input as:

res <- testDS_limma(d_counts, d_medians, design, contrast, block_id, trend, weights, markers_to_test, min_cells, min_samples, plot, path)

Looking further into testDS_limma() it seems like it only usable for paired data:

estimate correlation between paired samples

(note: paired designs only; >2 measures per sample not allowed)

But somehow runs on my data without error? I must be overlooking something…