Errors in getSnpMats, getSnpMats10X, and getCellAlleleCount: unable to find an inherited method for function ‘seqnames’ for signature ‘"DataFrame"’ about honeybadger

[JEFworks]

It looks like the supplied snps object input to your getSnpMats(snps, bamFilename, indexFilename) call is unable to pull out the seqnames slot. The seqnames slot provides the chromosome name that can then be used to get the relevant pileup from your bam files. It is possible that Granges class has been updated. Can you please check that you can run the following? load(system.file("ExAC", "ExAC_chr1.RData", package = "HoneyBADGER")) seqnames(snps) factor-Rle of length 6415 with 1 run Lengths: 6415 Values : 1 Levels(85): 1 2 3 4 5 ... GL000193.1 GL000194.1 GL000225.1 GL000192.1 NC_007605 On Apr 11, 2019, at 7:49 PM, koetjen <[email protected]<mailto:[email protected]>> wrote: I'm getting an unexpected error running through the prepare data tutorial on my 10X dataset. Is there any advice on where the source of the error might be and how to tackle it? results <- getSnpMats(snps, bamFilename, indexFilename) Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘seqnames’ for signature ‘"DataFrame"’ results <- getCellAlleleCount(snps, bamFilename, indexFilename, cellBarcodename) Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘seqnames’ for signature ‘"DataFrame"’ results <- getSnpMats10X(snps, bamFilename, indexFilename) Error in (function (classes, fdef, mtable) : unable to find an inherited method for function ‘seqnames’ for signature ‘"DataFrame"’ session_info() ─ Session info ─────────────────────────────────────────────────────────────────────────────────────────────────────────────────── setting value version R version 3.5.2 (2018-12-20) os macOS Mojave 10.14.3 system x86_64, darwin15.6.0 ui RStudio language (EN) collate en_US.UTF-8 ctype en_US.UTF-8 tz America/Chicago date 2019-04-11 ─ Packages ─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────── package * version date lib source AnnotationDbi 1.44.0 2018-10-30 [1] Bioconductor ape 5.3 2019-03-17 [1] CRAN (R 3.5.2) assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.5.2) backports 1.1.3 2018-12-14 [1] CRAN (R 3.5.0) bibtex 0.4.2 2017-06-30 [1] CRAN (R 3.5.0) Biobase * 2.42.0 2018-10-30 [1] Bioconductor BiocGenerics * 0.28.0 2018-10-30 [1] Bioconductor BiocManager 1.30.4 2018-11-13 [1] CRAN (R 3.5.0) BiocParallel * 1.16.6 2019-02-17 [1] Bioconductor biomaRt 2.38.0 2018-10-30 [1] Bioconductor Biostrings * 2.50.2 2019-01-03 [1] Bioconductor bit 1.1-14 2018-05-29 [1] CRAN (R 3.5.0) bit64 0.9-7 2017-05-08 [1] CRAN (R 3.5.0) bitops 1.0-6 2013-08-17 [1] CRAN (R 3.5.0) blob 1.1.1 2018-03-25 [1] CRAN (R 3.5.0) BSgenome 1.50.0 2018-10-30 [1] Bioconductor callr 3.2.0 2019-03-15 [1] CRAN (R 3.5.2) caTools 1.17.1.2 2019-03-06 [1] CRAN (R 3.5.2) cli 1.1.0 2019-03-19 [1] CRAN (R 3.5.2) cluster 2.0.7-1 2018-04-13 [1] CRAN (R 3.5.0) coda 0.19-2 2018-10-08 [1] CRAN (R 3.5.0) codetools 0.2-16 2018-12-24 [1] CRAN (R 3.5.2) colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.5.2) cowplot * 0.9.4 2019-01-08 [1] CRAN (R 3.5.2) crayon 1.3.4 2017-09-16 [1] CRAN (R 3.5.0) data.table 1.12.0 2019-01-13 [1] CRAN (R 3.5.2) DBI 1.0.0 2018-05-02 [1] CRAN (R 3.5.0) DelayedArray * 0.8.0 2018-10-30 [1] Bioconductor desc 1.2.0 2018-05-01 [1] CRAN (R 3.5.0) devtools * 2.0.1 2018-10-26 [1] CRAN (R 3.5.2) digest 0.6.18 2018-10-10 [1] CRAN (R 3.5.0) dplyr * 0.8.0.1 2019-02-15 [1] CRAN (R 3.5.2) fitdistrplus 1.0-14 2019-01-23 [1] CRAN (R 3.5.2) fs 1.2.7 2019-03-19 [1] CRAN (R 3.5.2) future 1.12.0 2019-03-08 [1] CRAN (R 3.5.2) future.apply 1.2.0 2019-03-07 [1] CRAN (R 3.5.2) gbRd 0.4-11 2012-10-01 [1] CRAN (R 3.5.0) gdata 2.18.0 2017-06-06 [1] CRAN (R 3.5.0) GenomeInfoDb * 1.18.2 2019-02-12 [1] Bioconductor GenomeInfoDbData 1.2.0 2019-03-12 [1] Bioconductor GenomicAlignments 1.18.1 2019-01-04 [1] Bioconductor GenomicFeatures 1.34.7 2019-03-25 [1] Bioconductor GenomicRanges * 1.34.0 2018-10-30 [1] Bioconductor ggplot2 * 3.1.1 2019-04-07 [1] CRAN (R 3.5.2) ggrepel 0.8.0 2018-05-09 [1] CRAN (R 3.5.0) ggridges 0.5.1 2018-09-27 [1] CRAN (R 3.5.0) globals 0.12.4 2018-10-11 [1] CRAN (R 3.5.0) glue 1.3.1 2019-03-12 [1] CRAN (R 3.5.2) gplots 3.0.1.1 2019-01-27 [1] CRAN (R 3.5.2) gridExtra 2.3 2017-09-09 [1] CRAN (R 3.5.0) gtable 0.3.0 2019-03-25 [1] CRAN (R 3.5.2) gtools 3.8.1 2018-06-26 [1] CRAN (R 3.5.0) hdf5r * 1.1.1 2019-03-26 [1] CRAN (R 3.5.2) hms 0.4.2 2018-03-10 [1] CRAN (R 3.5.0) HoneyBADGER * 0.1 2019-04-05 [1] Github (4eb9c06<https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_JEFworks_HoneyBADGER_commit_4eb9c06&d=DwMFaQ&c=WO-RGvefibhHBZq3fL85hQ&r=2gb0vmLv11Vi98WTAqlCXyDkhi11d9lKeGWDXEU-qNw&m=rOxdLDdHNQObTdYykJIR6BHr5F0tAkOBDtKrBmwhwOA&s=AH5BDRZoWBhcLwk-YusHIx3ipg89uo5VAW6kqkcjrgg&e=>) htmltools 0.3.6 2017-04-28 [1] CRAN (R 3.5.0) htmlwidgets 1.3 2018-09-30 [1] CRAN (R 3.5.0) httr 1.4.0 2018-12-11 [1] CRAN (R 3.5.0) ica 1.0-2 2018-05-24 [1] CRAN (R 3.5.0) igraph 1.2.4 2019-02-13 [1] CRAN (R 3.5.2) IRanges * 2.16.0 2018-10-30 [1] Bioconductor irlba 2.3.3 2019-02-05 [1] CRAN (R 3.5.2) jsonlite 1.6 2018-12-07 [1] CRAN (R 3.5.0) KernSmooth 2.23-15 2015-06-29 [1] CRAN (R 3.5.2) labeling 0.3 2014-08-23 [1] CRAN (R 3.5.0) lattice 0.20-38 2018-11-04 [1] CRAN (R 3.5.2) lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.5.2) listenv 0.7.0 2018-01-21 [1] CRAN (R 3.5.0) lmtest 0.9-36 2018-04-04 [1] CRAN (R 3.5.0) lsei 1.2-0 2017-10-23 [1] CRAN (R 3.5.0) magrittr 1.5 2014-11-22 [1] CRAN (R 3.5.0) MASS 7.3-51.4 2019-03-31 [1] CRAN (R 3.5.2) Matrix * 1.2-17 2019-03-22 [1] CRAN (R 3.5.2) matrixStats * 0.54.0 2018-07-23 [1] CRAN (R 3.5.0) memoise 1.1.0 2017-04-21 [1] CRAN (R 3.5.0) metap 1.1 2019-02-06 [1] CRAN (R 3.5.2) munsell 0.5.0 2018-06-12 [1] CRAN (R 3.5.0) nlme 3.1-137 2018-04-07 [1] CRAN (R 3.5.2) npsurv 0.4-0 2017-10-14 [1] CRAN (R 3.5.0) pbapply 1.4-0 2019-02-05 [1] CRAN (R 3.5.2) pillar 1.3.1 2018-12-15 [1] CRAN (R 3.5.0) pkgbuild 1.0.3 2019-03-20 [1] CRAN (R 3.5.2) pkgconfig 2.0.2 2018-08-16 [1] CRAN (R 3.5.0) pkgload 1.0.2 2018-10-29 [1] CRAN (R 3.5.0) plotly 4.8.0 2018-07-20 [1] CRAN (R 3.5.2) plyr 1.8.4 2016-06-08 [1] CRAN (R 3.5.0) png 0.1-7 2013-12-03 [1] CRAN (R 3.5.0) prettyunits 1.0.2 2015-07-13 [1] CRAN (R 3.5.0) processx 3.3.0 2019-03-10 [1] CRAN (R 3.5.2) progress 1.2.0 2018-06-14 [1] CRAN (R 3.5.0) ps 1.3.0 2018-12-21 [1] CRAN (R 3.5.0) purrr 0.3.2 2019-03-15 [1] CRAN (R 3.5.2) R.methodsS3 1.7.1 2016-02-16 [1] CRAN (R 3.5.0) R.oo 1.22.0 2018-04-22 [1] CRAN (R 3.5.0) R.utils 2.8.0 2019-02-14 [1] CRAN (R 3.5.2) R6 2.4.0 2019-02-14 [1] CRAN (R 3.5.2) RANN 2.6.1 2019-01-08 [1] CRAN (R 3.5.2) RColorBrewer 1.1-2 2014-12-07 [1] CRAN (R 3.5.0) Rcpp 1.0.1 2019-03-17 [1] CRAN (R 3.5.2) RCurl 1.95-4.12 2019-03-04 [1] CRAN (R 3.5.2) Rdpack 0.10-1 2018-10-04 [1] CRAN (R 3.5.0) remotes 2.0.2 2018-10-30 [1] CRAN (R 3.5.0) reshape2 * 1.4.3 2017-12-11 [1] CRAN (R 3.5.0) reticulate 1.11.1 2019-03-06 [1] CRAN (R 3.5.2) rjags 4-8 2018-10-19 [1] CRAN (R 3.5.0) rlang 0.3.4 2019-04-07 [1] CRAN (R 3.5.2) ROCR 1.0-7 2015-03-26 [1] CRAN (R 3.5.0) rprojroot 1.3-2 2018-01-03 [1] CRAN (R 3.5.0) Rsamtools * 1.34.1 2019-01-31 [1] Bioconductor RSQLite 2.1.1 2018-05-06 [1] CRAN (R 3.5.0) rstudioapi 0.10 2019-03-19 [1] CRAN (R 3.5.2) rsvd 1.0.0 2018-11-06 [1] CRAN (R 3.5.0) rtracklayer 1.42.2 2019-03-01 [1] Bioconductor Rtsne 0.15 2018-11-10 [1] CRAN (R 3.5.0) S4Vectors * 0.20.1 2018-11-09 [1] Bioconductor scales 1.0.0 2018-08-09 [1] CRAN (R 3.5.0) sctransform * 0.0.0.900 2019-04-08 [1] Github (satijalab/sctransform@4ab830c<https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_ChristophH_sctransform_commit_4ab830c&d=DwMFaQ&c=WO-RGvefibhHBZq3fL85hQ&r=2gb0vmLv11Vi98WTAqlCXyDkhi11d9lKeGWDXEU-qNw&m=rOxdLDdHNQObTdYykJIR6BHr5F0tAkOBDtKrBmwhwOA&s=VRaC0Z5-YneVy94rgPs0RVmYMl6Txf37-KhCP3lorT4&e=>) SDMTools 1.1-221 2014-08-05 [1] CRAN (R 3.5.0) sessioninfo 1.1.1 2018-11-05 [1] CRAN (R 3.5.0) Seurat * 3.0.0.9000 2019-03-30 [1] Github (satijalab/seurat@6d361c0<https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_satijalab_seurat_commit_6d361c0&d=DwMFaQ&c=WO-RGvefibhHBZq3fL85hQ&r=2gb0vmLv11Vi98WTAqlCXyDkhi11d9lKeGWDXEU-qNw&m=rOxdLDdHNQObTdYykJIR6BHr5F0tAkOBDtKrBmwhwOA&s=9nDfMCJo7RcSj_Fhn_jBoEBomiZt_ugN_XHdbdeB1ZA&e=>) stringi 1.4.3 2019-03-12 [1] CRAN (R 3.5.2) stringr 1.4.0 2019-02-10 [1] CRAN (R 3.5.2) SummarizedExperiment * 1.12.0 2018-10-30 [1] Bioconductor survival 2.44-1.1 2019-04-01 [1] CRAN (R 3.5.2) tibble 2.1.1 2019-03-16 [1] CRAN (R 3.5.2) tidyr 0.8.3 2019-03-01 [1] CRAN (R 3.5.2) tidyselect 0.2.5 2018-10-11 [1] CRAN (R 3.5.0) tsne 0.1-3 2016-07-15 [1] CRAN (R 3.5.0) usethis * 1.4.0 2018-08-14 [1] CRAN (R 3.5.0) VariantAnnotation * 1.28.13 2019-03-19 [1] Bioconductor viridisLite 0.3.0 2018-02-01 [1] CRAN (R 3.5.0) withr 2.1.2 2018-03-15 [1] CRAN (R 3.5.0) XML 3.98-1.19 2019-03-06 [1] CRAN (R 3.5.2) XVector * 0.22.0 2018-10-30 [1] Bioconductor yaml 2.2.0 2018-07-25 [1] CRAN (R 3.5.0) zlibbioc 1.28.0 2018-10-30 [1] Bioconductor zoo 1.8-5 2019-03-21 [1] CRAN (R 3.5.2) [1] /Library/Frameworks/R.framework/Versions/3.5/Resources/library — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub<https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_JEFworks_HoneyBADGER_issues_23&d=DwMFaQ&c=WO-RGvefibhHBZq3fL85hQ&r=2gb0vmLv11Vi98WTAqlCXyDkhi11d9lKeGWDXEU-qNw&m=rOxdLDdHNQObTdYykJIR6BHr5F0tAkOBDtKrBmwhwOA&s=ctkXuHVshvr9BPlFjuKXjj_0gEnbnAGcjUSFBSecsqY&e=>, or mute the thread<https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_notifications_unsubscribe-2Dauth_AIj2SApQ8sYBtd9dxei8MO9h20IzZ-5FPzks5vf8oJgaJpZM4crDu-2D&d=DwMFaQ&c=WO-RGvefibhHBZq3fL85hQ&r=2gb0vmLv11Vi98WTAqlCXyDkhi11d9lKeGWDXEU-qNw&m=rOxdLDdHNQObTdYykJIR6BHr5F0tAkOBDtKrBmwhwOA&s=NM_KPGFIt8eYs07hh-B63gV1PHQEh89K_uPzjP_H2YY&e=>.

[koetjen]

Thanks for the quick response. I can successfully run that. hmm, I wonder if that is related to this error:

testRanges <- GRanges("7", IRanges(start=1e5, width=1e6)) param <- ScanVcfParam(which=testRanges) vcf <- readVcf(vcfFiles[[i]], ref.genome, param=param)

Error: scanVcf: scanVcf: scanTabix: '7' not present in tabix index path: /Users/karolynoetjen/Link_Lab/3_scRNAseq/MDS/germline_VCF/H_QH-769516-1423296_2.annotated.coding_variant_filtered.vcf.gz index: /Users/karolynoetjen/Link_Lab/3_scRNAseq/MDS/germline_VCF/H_QH-769516-1423296_2.annotated.coding_variant_filtered.vcf.gz.tbi path: /Users/karolynoetjen/Link_Lab/3_scRNAseq/MDS/germline_VCF/H_QH-769516-1423296_2.annotated.coding_variant_filtered.vcf.gz

testRanges <- GRanges("chr7", IRanges(start=1e5, width=1e6)) param <- ScanVcfParam(which=testRanges) vcf <- readVcf(vcfFiles[[i]], ref.genome, param=param)

This is rest of my snps script and str: snps <- rowData(vcf) info <- info(vcf) maf <- info[, 'AF'] vi <- sapply(maf, function(x) any(x > maft)) snps <- snps[vi,] vi <- width(snps$REF) == 1 vi <- snps@listData$ALT@unlistData@ranges@width == 1 snps <- snps[vi,] str(snps) Formal class 'DataFrame' [package "S4Vectors"] with 6 slots ..@ rownames : chr [1:20] "chr7:246436_A/G" "chr7:291120_G/A" "chr7:291461_C/A" "chr7:550562_A/G" ... ..@ nrows : int 20 ..@ listData :List of 5 .. ..$ paramRangeID: Factor w/ 0 levels: NA NA NA NA NA NA NA NA NA NA ... .. ..$ REF :Formal class 'DNAStringSet' [package "Biostrings"] with 5 slots .. .. .. ..@ pool :Formal class 'SharedRaw_Pool' [package "XVector"] with 2 slots .. .. .. .. .. ..@ xp_list :List of 1 .. .. .. .. .. .. ..$ :<externalptr> .. .. .. .. .. ..@ .link_to_cached_object_list:List of 1 .. .. .. .. .. .. ..$ :<environment: 0x112b939d8> .. .. .. ..@ ranges :Formal class 'GroupedIRanges' [package "XVector"] with 7 slots .. .. .. .. .. ..@ group : int [1:20] 1 1 1 1 1 1 1 1 1 1 ... .. .. .. .. .. ..@ start : int [1:20] 4 2 1 4 4 1 2 4 4 2 ... .. .. .. .. .. ..@ width : int [1:20] 1 1 1 1 1 1 1 1 1 1 ... .. .. .. .. .. ..@ NAMES : NULL .. .. .. .. .. ..@ elementType : chr "ANY" .. .. .. .. .. ..@ elementMetadata: NULL .. .. .. .. .. ..@ metadata : list() .. .. .. ..@ elementType : chr "DNAString" .. .. .. ..@ elementMetadata: NULL .. .. .. ..@ metadata : list() .. ..$ ALT :Formal class 'DNAStringSetList' [package "Biostrings"] with 5 slots .. .. .. ..@ unlistData :Formal class 'DNAStringSet' [package "Biostrings"] with 5 slots .. .. .. .. .. ..@ pool :Formal class 'SharedRaw_Pool' [package "XVector"] with 2 slots .. .. .. .. .. .. .. ..@ xp_list :List of 1 .. .. .. .. .. .. .. .. ..$ :<externalptr> .. .. .. .. .. .. .. ..@ .link_to_cached_object_list:List of 1 .. .. .. .. .. .. .. .. ..$ :<environment: 0x112b939d8> .. .. .. .. .. ..@ ranges :Formal class 'GroupedIRanges' [package "XVector"] with 7 slots .. .. .. .. .. .. .. ..@ group : int [1:20] 1 1 1 1 1 1 1 1 1 1 ... .. .. .. .. .. .. .. ..@ start : int [1:20] 1 37 38 39 40 41 42 43 44 45 ... .. .. .. .. .. .. .. ..@ width : int [1:20] 1 1 1 1 1 1 1 1 1 1 ... .. .. .. .. .. .. .. ..@ NAMES : NULL .. .. .. .. .. .. .. ..@ elementType : chr "ANY" .. .. .. .. .. .. .. ..@ elementMetadata: NULL .. .. .. .. .. .. .. ..@ metadata : list() .. .. .. .. .. ..@ elementType : chr "DNAString" .. .. .. .. .. ..@ elementMetadata: NULL .. .. .. .. .. ..@ metadata : list() .. .. .. ..@ partitioning :Formal class 'PartitioningByEnd' [package "IRanges"] with 5 slots .. .. .. .. .. ..@ end : int [1:20] 1 2 3 4 5 6 7 8 9 10 ... .. .. .. .. .. ..@ NAMES : NULL .. .. .. .. .. ..@ elementType : chr "ANY" .. .. .. .. .. ..@ elementMetadata: NULL .. .. .. .. .. ..@ metadata : list() .. .. .. ..@ elementType : chr "DNAStringSet" .. .. .. ..@ elementMetadata: NULL .. .. .. ..@ metadata : list() .. ..$ QUAL : num [1:20] 477 552 691 373 237 ... .. ..$ FILTER : chr [1:20] "." "." "." "." ... ..@ elementType : chr "ANY" ..@ elementMetadata: NULL ..@ metadata : list()

[koetjen]

It looks like GenomicRanges deprecated seqnames in the most recent version

CHANGES IN VERSION 1.34.0

DEPRECATED AND DEFUNCT

o Deprecate several RangedData methods: seqinfo, seqinfo<-, seqnames, and
  findOverlaps#RangedData#GenomicRanges

[koetjen]

I apologize, closing this issue was inadvertent.

[JEFworks]

Thanks for the excellent debugging.

As you've seen, the getSnpMats, getCellAlleleCount, and getSnpMats10X functions are using the information provided in your snps object to obtain pileups in your provided bams with help from your provided index files.

As a result, the chromosome names in your snps object must match those used in your bams and subsequently your index files. By default, the chromosome names are provided in NCBI/Ensemble style (ie. without the 'chr'):

> library(HoneyBADGER)
> load(system.file("ExAC", "ExAC_chr1.RData", package = "HoneyBADGER"))
> seqlevelsStyle(snps)
[1] "NCBI"    "Ensembl"

If your snps object does not have 'chr' in its seqname but your reads are indexed with chr in their chromosome names (because you are using a different style), you will get errors as you've seen. Luckily, you can toggle between styles as follows:

> seqlevelsStyle(snps) <- "UCSC"
> seqlevels(snps)
 [1] "chr1"       "chr2"       "chr3"       "chr4"       "chr5"      
 [6] "chr6"       "chr7"       "chr8"       "chr9"       "chr10"     
[11] "chr11"      "chr12"      "chr13"      "chr14"      "chr15"     
[16] "chr16"      "chr17"      "chr18"      "chr19"      "chr20"     
[21] "chr21"      "chr22"      "chrX"       "chrY"       "chrM"      
[26] "GL000207.1" "GL000226.1" "GL000229.1" "GL000231.1" "GL000210.1"
[31] "GL000239.1" "GL000235.1" "GL000201.1" "GL000247.1" "GL000245.1"
[36] "GL000197.1" "GL000203.1" "GL000246.1" "GL000249.1" "GL000196.1"
[41] "GL000248.1" "GL000244.1" "GL000238.1" "GL000202.1" "GL000234.1"
[46] "GL000232.1" "GL000206.1" "GL000240.1" "GL000236.1" "GL000241.1"
[51] "GL000243.1" "GL000242.1" "GL000230.1" "GL000237.1" "GL000233.1"
[56] "GL000204.1" "GL000198.1" "GL000208.1" "GL000191.1" "GL000227.1"
[61] "GL000228.1" "GL000214.1" "GL000221.1" "GL000209.1" "GL000218.1"
[66] "GL000220.1" "GL000213.1" "GL000211.1" "GL000199.1" "GL000217.1"
[71] "GL000216.1" "GL000215.1" "GL000205.1" "GL000219.1" "GL000224.1"
[76] "GL000223.1" "GL000195.1" "GL000212.1" "GL000222.1" "GL000200.1"
[81] "GL000193.1" "GL000194.1" "GL000225.1" "GL000192.1" "NC_007605" 

Hope that helps.

• Jean

[koetjen]

I found that the problem arose from the initial structure of snps
In the tutorial:
snps <- rowData(vcf) # result is a data frame, with incorrect results for seqlevels() and does not contain seqnames()
Instead, it needs to be:
snps <- rowRanges(vcf) # result is a GRanges object, so seqlevels() matches the expected output and contains seqnames()
Thanks for all your help, Kary