6.2 Multi-group comparison

This section compares normpatch and TCC for a multi-group comparison. The example uses hypoData_mg, a simulated dataset distributed with TCC. It contains three groups with three replicates per group.

data(hypoData_mg)
group <- factor(c(1, 1, 1, 2, 2, 2, 3, 3, 3))
exp_design <- data.frame(group = group)

6.2.1 edgeR screening

For this multi-group example, both packages use matching TMM and edgeR settings during DEGES. Therefore, the estimated normalization factors are identical, as in the two-group TMM and edgeR comparison.

tcc_tmm <- new("TCC", hypoData_mg, group)
tcc_tmm <- calcNormFactors(
    tcc_tmm,
    norm.method = "tmm",
    test.method = "edger",
    iteration = 3,
    FDR = 0.1,
    floorPDEG = 0.05
)

tcc_tmm$norm.factors

##   G1_rep1   G1_rep2   G1_rep3   G2_rep1   G2_rep2   G2_rep3   G3_rep1   G3_rep2   G3_rep3 
## 1.0151652 0.9093168 0.7822821 0.8339844 1.1778906 1.2029681 1.0211199 1.0338776 1.0233953

x_tmm <- newSeqCountData(hypoData_mg, exp_design = exp_design)
x_tmm <- calc_nf(
    x_tmm,
    norm_func = norm_tmm,
    test_func = test_edger,
    iter = 3,
    p_cutoff = 0,
    q_cutoff = 0.1,
    min_rm_genes = 0.05
)

x_tmm@meta$nf 

##   G1_rep1   G1_rep2   G1_rep3   G2_rep1   G2_rep2   G2_rep3   G3_rep1   G3_rep2   G3_rep3 
## 1.0151652 0.9093168 0.7822821 0.8339844 1.1778906 1.2029681 1.0211199 1.0338776 1.0233953

tcc_tmm$norm.factors - x_tmm@meta$nf

## G1_rep1 G1_rep2 G1_rep3 G2_rep1 G2_rep2 G2_rep3 G3_rep1 G3_rep2 G3_rep3 
##       0       0       0       0       0       0       0       0       0

6.2.2 DESeq2 screening

For the multi-group comparison, the DESeq2-based DEGES workflow in TCC uses RLE normalization and an LRT to screen candidate DEGs. This matches the default normpatch DESeq2 workflow, so the normalization factors are expected to be identical with matching DEGES settings.

tcc_deseq2 <- new("TCC", hypoData_mg, group)
tcc_deseq2 <- calcNormFactors(
    tcc_deseq2,
    norm.method = "deseq2",
    test.method = "deseq2",
    iteration = 3,
    FDR = 0.1,
    floorPDEG = 0.05
)
tcc_deseq2$norm.factors

##   G1_rep1   G1_rep2   G1_rep3   G2_rep1   G2_rep2   G2_rep3   G3_rep1   G3_rep2   G3_rep3 
## 1.0115865 0.9061116 0.7945736 0.8290947 1.2077930 1.1774576 0.9939773 1.0175544 1.0618513

x_deseq2 <- newSeqCountData(hypoData_mg, exp_design = exp_design)
x_deseq2 <- calc_nf(
    x_deseq2,
    norm_func = norm_rle,
    test_func = test_deseq2,
    iter = 3,
    p_cutoff = 0,
    q_cutoff = 0.1,
    min_rm_genes = 0.05
)
x_deseq2@meta$nf

##   G1_rep1   G1_rep2   G1_rep3   G2_rep1   G2_rep2   G2_rep3   G3_rep1   G3_rep2   G3_rep3 
## 1.0115865 0.9061116 0.7945736 0.8290947 1.2077930 1.1774576 0.9939773 1.0175544 1.0618513

tcc_deseq2$norm.factors - x_deseq2@meta$nf

## G1_rep1 G1_rep2 G1_rep3 G2_rep1 G2_rep2 G2_rep3 G3_rep1 G3_rep2 G3_rep3 
##       0       0       0       0       0       0       0       0       0