Design model matrix with limma with cell different proportions estimates and batch effect

Question

[ij283]

Hello,

I am trying to run DEG while taking into account different life stages proportion for 6 samples in total, three replicates per conditions with 4 life stages (SA, SB, SC, SD) and batch effect (W_1). For that, I am using the functions below which seems only to deal with batch effect W_1. Can someone please help me to modify the code to be able to include the different estimates of life stages proportions before running the DE analysis ? All the codes works perfectly fine without A,B,C,D into the model matrix, hence I am asking for help. Many thanks in advance Tables I am trying to include in the model

---

A   B   C   D

1 0.163312321 0.106152566 0 0.005498893 2 2.409388389 1.040800283 0 0.106134021 3 0.01486834 0.010928011 0 0.000275666 4 1.57134067 0.276429515 0.110433492 0 5 2.186436606 0.417265851 0.172594234 0 6 1.706342807 0.339232907 0.113947311 0

Code for DE analysis starts here

batch.effect <- genes_batch$W `
       W_1

[1,] -0.6048296 [2,] 0.2495991 [3,] 0.3682783 [4,] 0.2955617 [5,] 0.2338163 [6,] -0.5424258

`design <- condition2design(condition = metatabl$label, batch.effect = batch.effect)

This is the condition2design function code I am using

condition2design <- function(condition,batch.effect=NULL){   design.data <- data.frame(condition=condition)  

if(!is.null(batch.effect)){ colnames(batch.effect) <- paste0('batch',1:ncol(batch.effect)) design.data <- data.frame(design.data,batch.effect) } design.fomula <- as.formula(paste0('~0+',paste0(colnames(design.data),collapse = '+'))) design <- model.matrix(design.fomula, data=design.data) idx <- grep('condition',colnames(design)) design.idx <- colnames(design)
design.idx <- substr(design.idx,nchar('condition')+1,nchar(design.idx))
colnames(design)[idx] <- design.idx[idx] design }

I then run the DE as follow:

if(DE_pipeline == 'limma'){

----->> limma pipeline result <- limma.pipeline(dge = genes_dge,

                             design = design,
                             deltaPS = NULL,
                             contrast = contrast,
                             diffAS = F,
                             adjust.method = pval_adj_method) }

if(DE_pipeline == 'glmQL'){

----->> edgeR glmQL pipeline result <- edgeR.pipeline(dge = genes_dge,

                             design = design,
                             deltaPS = NULL,
                             contrast = contrast,
                             diffAS = F,
                             method = 'glmQL',
                             adjust.method = pval_adj_method) }

if(DE_pipeline == 'glm'){

----->> edgeR glm pipeline result <- edgeR.pipeline(dge = genes_dge,

                             design = design,
                             deltaPS = NULL,
                             contrast = contrast,
                             diffAS = F,
                             method = 'glm',
                             adjust.method = pval_adj_method) }

please help with code amendments

[Gordon Smyth]

Cross-posted https://support.bioconductor.org/p/9158673/