[BioC] it is strange that I cannot get the fpkm by using easyRNASeq

[Fuyan [guest]]

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.