Bioconductor / Rsamtools

Binary alignment (BAM), FASTA, variant call (BCF), and tabix file import
https://bioconductor.org/packages/Rsamtools
Other
27 stars 27 forks source link

indexFa breaks FaFile object when using getSeq #5

Closed FelixErnst closed 5 years ago

FelixErnst commented 5 years ago

Hi,

I have noticed some odd behavior using the getSeq function with a FaFile object as input.

fa <- FaFile("test.fasta")
fa
> class: FaFile 
> path: test.fasta
> index: test.fasta.fai
> isOpen: FALSE 
> yieldSize: NA 
getSeq(FaFile(fasta),rtracklayer::import.gff3(gff))
> A DNAStringSet instance of length 1
>      width seq                          names               
> [1]     5 AGCTA                         chr1

The index file is automatically created, when I run getSeq.

However, If I now run the indexFa function, the FaFile object becomes unusable

indexFa(fa)
> class: FaFile 
> path: test.fasta
> index: test.fasta
> isOpen: FALSE 
> yieldSize: NA 
getSeq(fa,GRanges("chr1",IRanges(1,5)))
> Error in value[[3L]](cond) :  record 1 (chr1:1-5) failed
>   file: test.fasta

I noticed that the pointer for the index does revert to the fasta file itself, which might be the underlying problem. test.zip

I probably don't need to run the indexFa function to retrieve sequences. However, I noticed this behavior and just wanted to let someone know.

sessionInfo()
R version 3.5.2 (2018-12-20)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows >= 8 x64 (build 9200)

Matrix products: default

locale:
[1] LC_COLLATE=German_Germany.1252  LC_CTYPE=German_Germany.1252    LC_MONETARY=German_Germany.1252 LC_NUMERIC=C                   
[5] LC_TIME=German_Germany.1252    

attached base packages:
[1] stats4    parallel  stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] RNAmodR_0.0.2        Rsamtools_1.34.0     Modstrings_0.99.12   Biostrings_2.50.1    XVector_0.22.0       GenomicRanges_1.34.0
 [7] GenomeInfoDb_1.18.1  IRanges_2.16.0       S4Vectors_0.20.1     BiocGenerics_0.28.0 

loaded via a namespace (and not attached):
  [1] assertive.base_0.0-7        colorspace_1.3-2            colorRamps_2.3              biovizBase_1.30.1          
  [5] htmlTable_1.13.1            base64enc_0.1-3             dichromat_2.0-0             rstudioapi_0.8             
  [9] assertive.sets_0.0-3        bit64_0.9-7                 AnnotationDbi_1.44.0        assertive.data.uk_0.0-2    
 [13] codetools_0.2-15            splines_3.5.2               knitr_1.21                  Formula_1.2-3              
 [17] assertive_0.3-5             assertive.data.us_0.0-2     cluster_2.0.7-1             compiler_3.5.2             
 [21] httr_1.4.0                  backports_1.1.3             assertthat_0.2.0            Matrix_1.2-15              
 [25] lazyeval_0.2.1              acepack_1.4.1               htmltools_0.3.6             prettyunits_1.0.2          
 [29] tools_3.5.2                 bindrcpp_0.2.2              gtable_0.2.0                glue_1.3.0                 
 [33] GenomeInfoDbData_1.2.0      reshape2_1.4.3              dplyr_0.7.8                 Rcpp_1.0.0                 
 [37] Biobase_2.42.0              gdata_2.18.0                rtracklayer_1.42.1          assertive.files_0.0-2      
 [41] assertive.datetimes_0.0-2   assertive.models_0.0-2      xfun_0.4                    stringr_1.3.1              
 [45] ensembldb_2.6.3             gtools_3.8.1                XML_3.98-1.16               zlibbioc_1.28.0            
 [49] scales_1.0.0                BSgenome_1.50.0             VariantAnnotation_1.28.7    hms_0.4.2                  
 [53] ProtGenerics_1.14.0         SummarizedExperiment_1.12.0 AnnotationFilter_1.6.0      RColorBrewer_1.1-2         
 [57] assertive.matrices_0.0-2    assertive.strings_0.0-3     yaml_2.2.0                  curl_3.2                   
 [61] memoise_1.1.0               gridExtra_2.3               ggplot2_3.1.0               biomaRt_2.38.0             
 [65] rpart_4.1-13                latticeExtra_0.6-28         stringi_1.2.4               RSQLite_2.1.1              
 [69] checkmate_1.8.5             GenomicFeatures_1.34.1      caTools_1.17.1.1            BiocParallel_1.16.5        
 [73] rlang_0.3.0.1               pkgconfig_2.0.2             matrixStats_0.54.0          bitops_1.0-6               
 [77] lattice_0.20-38             assertive.data_0.0-3        ROCR_1.0-7                  purrr_0.2.5                
 [81] bindr_0.1.1                 GenomicAlignments_1.18.1    htmlwidgets_1.3             assertive.properties_0.0-4 
 [85] bit_1.1-14                  tidyselect_0.2.5            assertive.code_0.0-3        plyr_1.8.4                 
 [89] magrittr_1.5                R6_2.3.0                    gplots_3.0.1                Hmisc_4.1-1                
 [93] DelayedArray_0.8.0          DBI_1.0.0                   pillar_1.3.1                foreign_0.8-71             
 [97] assertive.numbers_0.0-2     survival_2.43-3             RCurl_1.95-4.11             nnet_7.3-12                
[101] tibble_1.4.2                crayon_1.3.4                assertive.types_0.0-3       KernSmooth_2.23-15         
[105] progress_1.2.0              grid_3.5.2                  data.table_1.11.8           blob_1.1.1                 
[109] digest_0.6.18               munsell_0.5.0               assertive.reflection_0.0-4  Gviz_1.26.4