con vs treated group miRNA comparison using limma
1
0
Entering edit mode
@d59ac2ea
Last seen 16 months ago
United States

Hi, I am using limma and EdgR to compare two groups. not sure which test I should use? Which test should be used to compare to variables? when we could use fdr, treat, glmQLFT or any other for unpaired comparison? Please correct me if I am using the wrong coding.

type <-factor(pheno_frame$Name)
design <- model.matrix(~ 0 + type, data=pheno_frame)
# Filtering, normalization, and transformation
keep <- rowSums(cpm(counts)>1.19) >= 8    #filter cpm>1 (50% of the samples)
x <- DGEList(counts = counts)
x <- x[keep, , keep.lib.sizes=FALSE]
dim(x)
cpmCutoff <- round(1e7/mean(colSums(counts)), 2)
grpNum <- 2
x <- x[rowSums(edgeR::cpm(x) > cpmCutoff) >= grpNum, ]
dim(x)
x <- calcNormFactors(x)
range(x$samples$norm.factors)
logcpm <- edgeR::cpm(x, log = TRUE, prior.count = 0.2)
boxplot(logcpm, las = 2)
x <- calcNormFactors(x,method="TMM") 
range(x$samples$norm.factors)
logcpm <- edgeR::cpm(x, log = TRUE, prior.count = 0.2)
boxplot(logcpm, las = 2)
# Differential miRNA between Test and Con
contrasts <-makeContrasts(Test_Con=typeTest-typecon,levels=design)
type <- c("Test","con")[as.factor(pheno_frame$Name)]
combo <- interaction(pheno_frame$HTG.EdgeSeq.Run.ID, type)
samp = pheno_frame$HTG.EdgeSeq.Run.ID
design <- model.matrix(~0+type )
pheno_frame$Name <- paste(pheno_frame$Name,sep = "_")
v <- voomWithQualityWeights(counts,plot=TRUE,normalize="none", design=design)
cpm <- 2^(v$E)
fit <- lmFit(v,design)
fit2 <-contrasts.fit(fit,contrasts)
fit2 <- eBayes(fit2)
Results <-topTable(fit2, coef="Test_Con",number=Inf, adjust="fdr",sort.by="P")
par(mfrow=c(1,2))
#v <- voom(d, design, plot=TRUE)
#v
array.weights <- arrayWeights(v, design)
fitw <- lmFit(v, design)#, weights = array.weights)
fit2 <- contrasts.fit(fitw, contrasts)
fit2 <- eBayes(fit2)
vfit <- lmFit(v, design)
vfit <- contrasts.fit(vfit, contrasts)
efit <- eBayes(vfit)
plotSA(fit2, main="Final model: Mean-variance trend")
#summary(decideTests(fit2))
#tfit <- treat(fit2, lfc=1)
#dt <- decideTests(tfit)
#summary(dt)
Results <- topTreat(fit2, coef="Test_Con", n=Inf, adjust="fdr",sort.by="P")
mirTarRnaSeq miRNAData • 608 views
ADD COMMENT
0
Entering edit mode
@gordon-smyth
Last seen 2 hours ago
WEHI, Melbourne, Australia

limma and edgeR provide a number of different analysis pipelines. The pipelines all perform well so you should just choose one instead of trying to mix them.

The packages provide a number of different analysis options, but you don't need to apply every option in every analysis. For any particular analysis, you only use the options necessary for your data. Your code is currently far too complicated and repetitive with lots of unnecessary complications. Just follow one of the documented examples and keep the analysis as simple as possible. Comparing two groups is actually very easy and can be done in a few lines of code.

If you want to apply the limma-voom pipeline, then a very complete worked example is provided here: limma-voom workflow.

ADD COMMENT

Login before adding your answer.

Traffic: 512 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6