reutrn "na" values for basemeanA and basemeanB with DESeq2 run

[young]

hello,

I am doing DE analysis with DESeq2 method and it is perfectly working fine.

But when i checked its results, it return "na" for basemeanA and basemeanB.

when i tried with previous Trinity,Deseq2 version,  basemeanA and basemeanB have values for each group,

but with recent version, no values return. please help with this?  

Thank you,

#version info##########################################################################

Trinity v2.4.0

edgeR : v3.18.1

DESeq2 : v1.16.1

#script run###########################################################################

perl trinityrnaseq-Trinity-v2.4.0/Analysis/DifferentialExpression/run_DE_analysis.pl           \
                                                               --matrix trans_counts.counts.matrix                       \
                                                               --method DESeq2                                                  \
                                                               --samples_file sample.txt                                       \
                                                               --min_reps_min_cpm 2,1                                      \
                                                               --output out_results                                               \
                                                               --contrasts contrast.txt                                          \
                                                            

#Rscript#############################################################################

data = read.table("/rna/SEVERANCE_YangWooKyeom_ampliseq/2017-08-09_human_66sp_deg_2/analysis/trans_counts.counts.matrix", header=T, row.names=1, com='')
col_ordering = c(1,2,3,5,6,19,20,21,22,8,9,10,45,47,50,52,59,64,66)
rnaseqMatrix = data[,col_ordering]
rnaseqMatrix = round(rnaseqMatrix)
rnaseqMatrix = rnaseqMatrix[rowSums(cpm(rnaseqMatrix) > 1) >= 2,]
conditions = data.frame(conditions=factor(c(rep("GroupA", 9), rep("GroupC", 10))))
rownames(conditions) = colnames(rnaseqMatrix)
ddsFullCountTable <- DESeqDataSetFromMatrix(
    countData = rnaseqMatrix,
    colData = conditions,
    design = ~ conditions)
dds = DESeq(ddsFullCountTable)
contrast=c("conditions","GroupA","GroupC")
res = results(dds, contrast)
baseMeanA <- rowMeans(counts(dds, normalized=TRUE)[,colData(dds)$condition == "GroupA"])
baseMeanB <- rowMeans(counts(dds, normalized=TRUE)[,colData(dds)$condition == "GroupC"])
res = cbind(baseMeanA, baseMeanB, as.data.frame(res))
res = cbind(sampleA="GroupA", sampleB="GroupC", as.data.frame(res))
res$padj[is.na(res$padj)]  <- 1
write.table(as.data.frame(res[order(res$pvalue),]), file='trans_counts.counts.matrix.GroupA_vs_GroupC.DESeq2.DE_results', sep=' ', quote=FALSE)

#DEG output#########################################################################

sampleA sampleB baseMeanA baseMeanB baseMean log2FoldChange lfcSE stat pvalue padj
RPS26 GroupA GroupC NA NA 3243.52630878494 2.32009769400642 0.251327827261281 9.23136016925987 2.67216778069705e-20 4.51195529770696e-16
USP17L6P GroupA GroupC NA NA 643.911436330647 3.17671797847387 0.477328494423506 6.65520289609057 2.8290949916125e-11 1.95889334051027e-07
PPP1R9B GroupA GroupC NA NA 524.160967263814 -1.48204932121481 0.223716940534865 -6.624662923029 3.48041458189565e-11 1.95889334051027e-07
ZNF252P-AS1 GroupA GroupC NA NA 61.3339041685208 2.44410550299089 0.385029645695761 6.3478371868595 2.18362996616062e-10 7.7375158489103e-07
ZC3HAV1 GroupA GroupC NA NA 305.866148167732 -0.93015669909755 0.147011043404723 -6.32712126623589 2.49777061169271e-10 7.7375158489103e-07
ZNF625 GroupA GroupC NA NA 321.279231278962 1.14775812573997 0.181829318415566 6.31228305611748 2.74948742040046e-10 7.7375158489103e-07
HTATSF1 GroupA GroupC NA NA 114.191844641627 -0.88251248819292 0.150607653152547 -5.85967890555365 4.63763013927075e-09 9.87032434253701e-06
ITPK1 GroupA GroupC NA NA 264.264741463003 -0.997008710470226 0.170322358592377 -5.85365725739105 4.80879125161566e-09 9.87032434253701e-06

[Brian Haas]

Hi,

If you can package up your inputs and send them to me, I'll take a look.

best,

~brian

--
You received this message because you are subscribed to the Google Groups "trinityrnaseq-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to trinityrnaseq-users+unsub...@googlegroups.com.
To post to this group, send email to trinityrnaseq-users@googlegroups.com.
Visit this group at https://groups.google.com/group/trinityrnaseq-users.
For more options, visit https://groups.google.com/d/optout.

--

--
Brian J. Haas
The Broad Institute
http://broadinstitute.org/~bhaas

 

[KARAN PATEL]

Same problem for me. So, what can i do to resolve it?

To post to this group, send email to trinityrn...@googlegroups.com.

Visit this group at https://groups.google.com/group/trinityrnaseq-users.
For more options, visit https://groups.google.com/d/optout.

[Brian Haas]

Hi,

NA values were usually tied to the issue of quantification having been done on different target assemblies (as opposed to a single assembly based on all combined reads).

Could that be the case here?

best,

~brian

To unsubscribe from this group and stop receiving emails from it, send an email to trinityrnaseq-u...@googlegroups.com.

To post to this group, send email to trinityrn...@googlegroups.com.
Visit this group at https://groups.google.com/group/trinityrnaseq-users.
For more options, visit https://groups.google.com/d/optout.

--

--
Brian J. Haas
The Broad Institute
http://broadinstitute.org/~bhaas

 

--

You received this message because you are subscribed to the Google Groups "trinityrnaseq-users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to trinityrnaseq-u...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/trinityrnaseq-users/70397974-c944-4faf-a413-1cd43783e2ee%40googlegroups.com.

[JM]

Hi Brian,

I am facing the same issue. I did have to assemble the target transcriptome including an extra set of reads that were not aligned - they were not part of the differential expression experiment, instead they were added to improve the target transcriptome assembly. Will having these NAs lead to meaningful consequences for downstream analyses or in interpreting results?

Thanks for your help.

Cheers,

Jasmine

[Brian Haas]

possibly.  There really shouldn't be NA values - it's a symptom of something not being done according to the protocol.

To view this discussion on the web visit https://groups.google.com/d/msgid/trinityrnaseq-users/326e36e1-977c-45f6-aa00-0abaa2583f8fn%40googlegroups.com.