[BioC] it is strange that I cannot get the fpkm by using easyRNASeq
Nicolas Delhomme
nicolas.delhomme at plantphys.umu.se
Tue Apr 9 13:58:06 CEST 2013
Dear Fuyan,
No worries.
You seem to have put the finger on a bug, I'll check that out. A new version of R and Bioc is available and any change I would do would go into that new version (will be 1.6.1 then).
As we discussed offline, remember that easyRNASeq won't give you FPKM sensus-stricto. What you really get are RPKMs.
Nico
---------------------------------------------------------------
Nicolas Delhomme
Nathaniel Street Lab
Department of Plant Physiology
Umeå Plant Science Center
Tel: +46 90 786 7989
Email: nicolas.delhomme at plantphys.umu.se
SLU - Umeå universitet
Umeå S-901 87 Sweden
---------------------------------------------------------------
On Apr 9, 2013, at 11:13 AM, Fuyan [guest] wrote:
>
> Dear Nico,
>
> Sorry to disturb you again.
>
> I am trying to use easyRNASeq to calculate fpkm values.
>
>
>
> I used the sample from the vignette, but all I get are "NA".
>
> Below are the details:
>
>
>
>
>> count.table <- easyRNASeq(system.file(
> + "extdata",
> + package="RnaSeqTutorial"),
> + organism="Dmelanogaster",
> + readLength=30L,
> + annotationMethod="rda",
> + annotationFile=system.file(
> + "data",
> + "gAnnot.rda",
> + package="RnaSeqTutorial"),
> + count="exons",
> + filenames=c("ACACTG.bam", "ACTAGC.bam",
> + "ATGGCT.bam", "TTGCGA.bam"),
> + normalize=TRUE
> + )
> Checking arguments...
> Fetching annotations...
> Summarizing counts...
> Processing ACACTG.bam
> Processing ACTAGC.bam
> Processing ATGGCT.bam
> Processing TTGCGA.bam
> Preparing output
> Normalizing counts
> Warning messages:
> 1: In easyRNASeq(system.file("extdata", package = "RnaSeqTutorial"), :
> There are 50573 features/exons defined in your annotation that overlap! This implies that some reads will be counted more than once! Is that really what you want?
> 2: In easyRNASeq(system.file("extdata", package = "RnaSeqTutorial"), :
> You enforce UCSC chromosome conventions, however the provided annotation is not compliant. Correcting it.
> 3: In fetchCoverage(rnaSeq, format = format, filename = filename, filter = filter, :
> You enforce UCSC chromosome conventions, however the provided alignments are not compliant. Correcting it.
> 4: In fetchCoverage(rnaSeq, format = format, filename = filename, filter = filter, :
> You enforce UCSC chromosome conventions, however the provided alignments are not compliant. Correcting it.
> 5: In fetchCoverage(rnaSeq, format = format, filename = filename, filter = filter, :
> You enforce UCSC chromosome conventions, however the provided alignments are not compliant. Correcting it.
> 6: In fetchCoverage(rnaSeq, format = format, filename = filename, filter = filter, :
> You enforce UCSC chromosome conventions, however the provided alignments are not compliant. Correcting it.
>
>
>
>
>> head(count.table)
> ACACTG.bam ACTAGC.bam ATGGCT.bam TTGCGA.bam
> CG11023:1 NA NA NA NA
> CG11023:2 NA NA NA NA
> CG11023:3 NA NA NA NA
> CG2671:1 NA NA NA NA
> CG2671:2 NA NA NA NA
> CG2671:3 NA NA NA NA
>
>
>
> while the easyRNASeq works well without normalize=TRUE:
>> head(count.table)
> ACACTG.bam ACTAGC.bam ATGGCT.bam TTGCGA.bam
> CG11023:1 0 0 0 0
> CG11023:2 0 0 0 0
> CG11023:3 0 0 0 1
> CG2671:1 0 0 0 0
> CG2671:2 1 0 0 1
> CG2671:3 13 8 11 12
>
>
>
>
>
>
> I do not know what is wrong with my program.
>
> When I am trying to use my data, I also only get "NA" in the result.
>
> -- output of sessionInfo():
>
>> sessionInfo()
> R version 2.15.1 (2012-06-22)
> Platform: i386-pc-mingw32/i386 (32-bit)
>
> locale:
> [1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 LC_MONETARY=English_United States.1252
> [4] LC_NUMERIC=C LC_TIME=English_United States.1252
>
> attached base packages:
> [1] parallel stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] BSgenome.Dmelanogaster.UCSC.dm3_1.3.19 RnaSeqTutorial_0.0.11 easyRNASeq_1.4.2
> [4] ShortRead_1.16.4 latticeExtra_0.6-24 RColorBrewer_1.0-5
> [7] Rsamtools_1.10.2 DESeq_1.10.1 lattice_0.20-14
> [10] locfit_1.5-8 BSgenome_1.26.1 GenomicRanges_1.10.7
> [13] Biostrings_2.26.3 IRanges_1.16.6 edgeR_3.0.8
> [16] limma_3.14.4 biomaRt_2.14.0 Biobase_2.18.0
> [19] genomeIntervals_1.14.0 BiocGenerics_0.4.0 intervals_0.14.0
>
> loaded via a namespace (and not attached):
> [1] annotate_1.36.0 AnnotationDbi_1.20.7 bitops_1.0-5 DBI_0.2-5 genefilter_1.40.0 geneplotter_1.36.0 grid_2.15.1
> [8] hwriter_1.3 RCurl_1.95-4.1 RSQLite_0.11.2 splines_2.15.1 stats4_2.15.1 survival_2.37-4 tools_2.15.1
> [15] XML_3.95-0.2 xtable_1.7-1 zlibbioc_1.4.0
>
> --
> Sent via the guest posting facility at bioconductor.org.
More information about the Bioconductor
mailing list