SEQanswers

Go Back   SEQanswers > Bioinformatics > Bioinformatics



Similar Threads
Thread Thread Starter Forum Replies Last Post
EDGE-pro into DESeq, DESeq run error? epistatic Bioinformatics 5 09-02-2015 05:32 AM
interpretation of fastQC result liaojinyue Bioinformatics 1 07-16-2014 01:30 PM
strange plotDispEsts result from Deseq sadiexiaoyu Bioinformatics 2 06-03-2013 01:25 PM
Strange Result in Deseq Output cburke11 Bioinformatics 2 12-04-2012 09:35 AM
interpreting DESeq result, padj - values Azazel Bioinformatics 0 10-06-2010 06:49 PM

Reply
 
Thread Tools
Old 08-04-2014, 05:14 PM   #1
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default DESeq result interpretation

Dear ALL,
this is my first time i use DESeq for two RNA-seq samples control and affected
for one gene using exon data
id baseMean baseMeanA baseMeanB foldChange log2FoldChange pval padj
1 122.3801454 86.64640789 158.113883 1.824817518 0.867752202 0.227324022 1
2 82.8516747 73.99729725 91.70605214 1.239316239 0.30954437 0.701941312 1
3 52.65192304 46.80170937 58.50213671 1.25 0.321928095 0.736690545 1
4 106.726871 98.03060747 115.4231346 1.177419355 0.235628248 0.754103212 1
5 34.94316814 49.3315315 20.55480479 0.416666667 -1.263034406 0.271418019 1
6 17.7087549 25.93067681 9.486832981 0.365853659 -1.450661409 0.375252652 1
7 28.46049894 28.46049894 28.46049894 1 3.20E-16 1 1
8 24.82387963 22.76839915 26.87936011 1.180555556 0.239465935 0.873379607 1
9 13.43968006 18.97366596 7.90569415 0.416666667 -1.263034406 0.521556208 1
10 48.3828482 49.3315315 47.4341649 0.961538462 -0.056583528 0.999156423 1
11 436.2362032 532.527558 339.9448485 0.638361045 -0.647555479 0.274798212 1


i need to know what type of info i can get from this result
i will be obliged if somebody will give a feedback. it will encourage me to explore my results for whole transcriptome data
huma Asif is offline   Reply With Quote
Old 08-05-2014, 02:09 AM   #2
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

It's unclear what you mean by "exon data" in this context, though I'm assuming that you mean "I used counts of the number of reads mapping to the exons as the counts for each gene".

The most important things to look at in such results are the adjusted p-values (the last column) and fold changes (the 5th or 6th column, depending on whether you prefer them log2 transformed). In your case, the snippet you posted contains no differentially expressed genes, at least when one uses a standard threshold for significance.

BTW, since you're starting out, you should use DESeq2, which has a number of improved features.
dpryan is offline   Reply With Quote
Old 08-05-2014, 05:59 AM   #3
velt
Member
 
Location: Paris

Join Date: Jun 2013
Posts: 10
Default

I advise you to sort your files in ascending order on the "adjusted p-value" column. Thus, the most significant genes are in the beginning of file. Generally a p-value is significant if it is less than 0.05. Then, as stated by dpryan, refer you to the "log2 fold-change" column. More the value is high (or low, less than 0), more the change of expression is important.
velt is offline   Reply With Quote
Old 08-05-2014, 06:02 AM   #4
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default thanks

yes you are right.
can you guide me where i can read more about interpretation.
i am not with statistical background so having problem in understanding in manual and paper they didnt write much.
huma Asif is offline   Reply With Quote
Old 08-05-2014, 06:03 AM   #5
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

in my case padjusted is all one and p value for any exon is not equal to 0.05
let me clear i am checking just one gene with 11exons
huma Asif is offline   Reply With Quote
Old 08-05-2014, 06:12 AM   #6
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

If you want to look at differential exon usage then you should use DEXseq rather than DESeq/DESeq2. If you want to look at differential expression, then each gene should only have one count associated with it (i.e., you need to sum the exons (depending on how you did the counting)).
dpryan is offline   Reply With Quote
Old 08-05-2014, 06:21 AM   #7
velt
Member
 
Location: Paris

Join Date: Jun 2013
Posts: 10
Default

If you want to learn more about DESeq or DESeq2, the documentation is very well done.

DESeq:

http://bioconductor.org/packages/rel.../doc/DESeq.pdf

DESeq2:

https://bioconductor.org/packages/de...doc/DESeq2.pdf

It is difficult with the little information that you give to explain the way to go. But the documentation is well, with the "pasilla dataset" that you can use in R to understand the different steps.
velt is offline   Reply With Quote
Old 08-05-2014, 06:32 AM   #8
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

ok so you mean DESeq work with gene and DEXseq with exon
if yes then i will install DEXseq and try again
i will also try DESeq by summing the counts for every exon into one gene
huma Asif is offline   Reply With Quote
Old 08-05-2014, 06:44 AM   #9
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Yes, that's exactly what I mean.
dpryan is offline   Reply With Quote
Old 08-05-2014, 06:25 PM   #10
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default Error

hi i tried to used DESeq with gene count and got this error after
countDataset<-estimateDispersions(countDataset,method='blind',sharingMode="fit-only")

tried this also but error persist
countDataset<-estimateDispersions(countDataset,method='blind')

Error in glm.fit(x = numeric(0), y = numeric(0), weights = NULL, start = c(0.1, :
object 'fit' not found
In addition: Warning messages:
1: In glm.fit(x = numeric(0), y = numeric(0), weights = NULL, start = c(0.1, :
no observations informative at iteration 1
2: glm.fit: algorithm did not converge
huma Asif is offline   Reply With Quote
Old 08-06-2014, 01:31 AM   #11
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Is there a reason you're using blind dispersion estimation? Aside from that, it looks like a row with 0 counts didn't get removed prior to testing.

BTW, as I mentioned before, I strongly recommend using DESeq2 instead.
dpryan is offline   Reply With Quote
Old 08-06-2014, 04:10 AM   #12
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

for DEXSeq i downloaded Homo_sapiens.gtf from UCSC but it terminate everytime on chr15
from ensembl it is downloaded successfully but now prob is chr in bam and 1,2,... chromosome notation
can you guide how to make it compatible
huma Asif is offline   Reply With Quote
Old 08-06-2014, 04:30 AM   #13
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Regarding the Ensembl annotation, you could just use sed or awk to change the chromosome names (note that the mitochondria are chrM in UCSC and MT in Ensembl, so just appending "chr" won't fix everything).

Regarding the DEXSeq error, it would be helpful to see the command issued and the exact error message.
dpryan is offline   Reply With Quote
Old 08-06-2014, 09:34 AM   #14
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

thanks i downloaded gtf from ucsc

well i am following instruction in
Inferring differential exon usage in RNA-Seq data with the DEXSeq package
done with
dexseq_prepare_annotation.py script and
python /path/to/library/DEXSeq/python_scripts/dexseq_count.py
Dmel_flattenend.gff untreated1.sam untreated1.counts

done with > sampleTable
countFile condition libType
A A.counts control single-end
B B.counts patient single-end


need help here
ecs <- read.HTSeqCounts(
+ sampleTable$countFile,
+ sampleTable,
+ "gene.gff" )

getting this error
Error: unexpected string constant in:
"sampleTable
"gene.gff""

Last edited by huma Asif; 08-06-2014 at 09:36 AM.
huma Asif is offline   Reply With Quote
Old 08-06-2014, 10:55 AM   #15
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

What version of DEXseq are you using? That entire method has been deprecated. Anyway, try:
Code:
flattenedfile="gene.gff"
instead of just "gene.gff".
dpryan is offline   Reply With Quote
Old 08-06-2014, 11:43 AM   #16
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

i am using
> packageVersion("DEXSeq")
[1] ‘1.10.8’
huma Asif is offline   Reply With Quote
Old 08-06-2014, 11:50 AM   #17
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

> ecs <- read.HTSeqCounts(
+ sampleTable$countFile,
+ sampleTable,
+ flattenedfile="gene.gff")
Error in checkAtAssignment("character", "annotationFile", "character") :
‘annotationFile’ is not a slot in class “character”
In addition: Warning message:
'newExonCountSet' is deprecated.
Use 'DEXSeqDataSet' instead.
See help("Deprecated")

Last edited by huma Asif; 08-06-2014 at 12:05 PM.
huma Asif is offline   Reply With Quote
Old 08-06-2014, 01:34 PM   #18
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

please help me in creating
DEXSeqDataSet
i am confused what shall i give in place of countFiles
where do i specify my files in this code



dxd = DEXSeqDataSetFromHTSeq(
countFiles,
sampleData=sampleTable,
design= sample + exon + condition:exon,
flattenedfile=flattenedFile )
huma Asif is offline   Reply With Quote
Old 08-07-2014, 01:50 AM   #19
dpryan
Devon Ryan
 
Location: Freiburg, Germany

Join Date: Jul 2011
Posts: 3,480
Default

Please read
Code:
help(DEXSeqDataSetFromHTSeq)
You'll see that the countFiles are the files containing counts produced by dexseq_count.py.
dpryan is offline   Reply With Quote
Old 08-07-2014, 05:27 AM   #20
huma Asif
Member
 
Location: Japan

Join Date: Oct 2010
Posts: 53
Default

thank you i read all that but do i need to write names of my count file or path where my files are my current diectory is where my gff and count files are
huma Asif is offline   Reply With Quote
Reply

Thread Tools

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off




All times are GMT -8. The time now is 09:39 PM.


Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2020, vBulletin Solutions, Inc.
Single Sign On provided by vBSSO