Dear Wolfgang, dear Simon,
Thanks for the explanation, I did not have understand well the purpose
of DEXSeq.
If you are interested, here are the cases in my dataset where genes
have more than 1 exon, but testGeneForDEU did not give a result:
*The simplest case:*
res1[1449:1451,]
geneID exonID dispersion_CR_est dispersion pvalue
EHI_031280:001 EHI_031280 E001 0.25177357 2.517736e-01 0.837492
EHI_031280:002 EHI_031280 E002 0.03580064 3.580064e-02 0.837492
EHI_031280:003 EHI_031280 E003 NA 1.000000e+08 NA
padjust
EHI_031280:001 0.9746982
EHI_031280:002 0.9746982
EHI_031280:003 NA
And the read counts per exon are:
exon_counts[1449:1451,]
exon_ID FKKM1_1 FPKM1_2 FPKM1_3 FPKM2_1 FPKM2_2 FPKM2_3
1449 EHI_031280:001 9 9 5 33 7 2
1450 EHI_031280:002 155 192 220 453 133 176
1451 EHI_031280:003 0 0 0 0 0 0
One exon in the gene has 0 for all FPKM => 10^8 dispersion and no
result for the test. I am not surprised, contrary to *this second
case*:
res1[1227:1228,]
geneID exonID dispersion_CR_est dispersion pvalue padjust
EHI_025270:001 EHI_025270 E001 NA 1.00000e+08 NA NA
EHI_025270:002 EHI_025270 E002 NA 2.03465e-02 NA NA
And the read counts per exon are:
exon_counts[1227:1228,]
exon_ID FKKM1_1 FPKM1_2 FPKM1_3 FPKM2_1 FPKM2_2 FPKM2_3
1227 EHI_025270:001 0 0 0 0 0 0
1228 EHI_025270:002 11900 3012 12136 10971 7387 5923
Why is there no result for the second exon?
Thanks,
Jane
"Dear Jane,
in the case of genes (or better term: loci) with only one exon, DEXSeq's
testGeneForDEU does not make sense - it is looking for exons within the
locus that behave differently from the average of all the others: there
are no others! Hence the parameters are unidentifiable, and the software
should not even look at them. We'll robustify the software for this case
in future versions (note that the package is still actively being
developed).
The case of the NA values that are not explained by there only being one
exon is more interesting. Can you send us the code to reproduce your
problem (this also requires the data, which you can send offline, and in
which you can for instance anonymize the gene and sample annotations).
Best wishes
Wolfgang
Hi Jane
On 08/08/2011 02:45 PM, Wolfgang Huber wrote:
>* in the case of genes (or better term: loci) with only one exon, DEXSeq's
*>* testGeneForDEU does not make sense - it is looking for exons within the
*>* locus that behave differently from the average of all the others: there
*>* are no others! Hence the parameters are unidentifiable, and the software
*>* should not even look at them.
*
Actually, it doesn't. This is why the p value is NA. Also note the
column 'testable' in the fData data frame: it is set to FALSE in such
cases to indicate that a test was not attempted. You will also find
testable=FALSE occasionally in genes with more than one exons, namely
whenever all exons or all exons but one have zero counts in all samples.
Cheers
Simon"
2011/7/29 Jane Merlevede <jane.merlev...@gmail.com>
> Hello,
>
> I am using DEXseq for differential analysis. I have posted some e-mails
> about it already on this list, but I have more questions !
> My dataset has 2 conditions with 3 biological replicates per condition:
>
> pData(ecs)
> sizeFactor condition replicate type
> Raman_1 0.9638822 nonvir 1 paired-end
> HM1_1 1.4000715 vir 1 paired-end
> Raman_2 1.0314049 nonvir 2 paired-end
> Raman_3 1.0734133 nonvir 3 paired-end
> HM1_2 0.7801269 vir 2 paired-end
> HM1_3 0.9488409 vir 3 paired-end
>
>
> I found the number of differentially expressed exons. I would like to know
> also which are down and which are up regulated. This information is given by
> foldChange, which is not provided in DEXSeq (res1).
>
> res1 <- DEUresultTable(ecs)
> head(res1)
> geneID exonID dispersion_CR_est dispersion pvalue
> EHI_000010:001 EHI_000010 E001 NA 0.02377712 NA
> EHI_000130:001 EHI_000130 E001 0.07466945 0.07466945 0.4867354
> EHI_000130:002 EHI_000130 E002 0.12897281 0.12897281 0.9427953
> EHI_000130:003 EHI_000130 E003 0.01960222 0.02031871 0.8432791
> EHI_000130:004 EHI_000130 E004 0.06733720 0.06733720 0.7042977
> EHI_000240:001 EHI_000240 E001 0.04375400 0.04375400 0.7493987
> padjust
> EHI_000010:001 NA
> EHI_000130:001 0.8639624
> EHI_000130:002 0.9927663
> EHI_000130:003 0.9756197
> EHI_000130:004 0.9503086
> EHI_000240:001 0.9611397
>
> So I would like to know, if it could be possible to have baseMean,
> baseMeanA, baseMeanB, foldChange, log2FoldCHange, resVarA and resVarB as in
> DESeq?
>
> Since this does not seem to be implemented in DEXSeq, I should start with
> the read counts per exon table to get it, but I can't find it neither. I
> think it is use in the estimateSizeFactors function, but I don't manage to
> access it. It's a pity, I will have to reconstruct the table...
>
>
> Thank you for your help,
> Jane Merlevède
>
[[alternative HTML version deleted]]
_______________________________________________
Bioc-sig-sequencing mailing list
Bioc-sig-sequencing@r-project.org
https://stat.ethz.ch/mailman/listinfo/bioc-sig-sequencing
