I'm conducting a DE analysis over normalized TPM, RPKM, or FPKM data. Raw counts are not available. The data is about the clinical outcome of a certain drug. I know this is not the ideal way to do this kind of analysis, but I have no choice. This complicates things for me, I'm not getting results and I just want to make sure my code is right. Here is my code, from one of the datasets, that has TPM data:
TPM = TPM[rowSums(TPM)>0,] thresh = TPM > 10 keep = rowSums(thresh) >= 12 table(keep) TPM = TPM[keep,] # Design and Contrast design = model.matrix(~ 0 + outcome, metadata) colnames(design)[c(1,2)] = c('NoResponse','Response') contrast = makeContrasts(NoResp_vs_Resp = NoResponse - Response, levels = colnames(design)) log_TPM = scale(log2(TPM + 0.1)) vfit = lmFit(log_TPM, design) vfit = contrasts.fit(vfit, contrasts = contrast) efit = eBayes(vfit, robust = TRUE, trend = TRUE) plotSA(efit, main = 'final model with TPM ~ Mean-Variance trend') summary(decideTests(efit)) DEG = topTable(efit, coef = 'NoResp_vs_Resp', p.value = 0.05, adjust.method = 'fdr', number = Inf)
I'm trying to get figures that are as close to the Bioconductors guide for limma.
Can anybody please just confirm wheather my code is alright, or is there any other way I could do this? regarding filtering and stuff like that, cause I cant use