Dear list, when I try to count human reads vs a human transcript database, I get the following error: library(GenomicFeatures) library(biomaRt) Here I construct my human transcript database: ensembl = useDataset("hsapiens_gene_ensembl",mart=useMart("ensembl")) res = getBM(attributes = "ensembl_transcript_id", filters = ("biotype","status"),values=list("protein_coding","known"),mart = ensembl) humanDb = makeTranscriptDbFromBiomart(biomart = "ensembl", dataset = "hsapiens_gene_ensembl",transcript_ids = as.character(res[,1])) tx = transcriptsBy(humanDb) The reads are an GappedAlignments object: reads GappedAlignments with 11791172 alignments and 0 elementMetadata cols: seqnames strand cigar qwidth start end <rle> <rle> <character> <integer> <integer> <integer> SRR015293.3 3 * 32M 32 186338939 186338970 SRR015293.5 16 * 32M 32 72094409 72094440 SRR015293.7 1 * 32M 32 159683461 159683492 SRR015293.9 4 * 32M 32 155529684 155529715 SRR015293.10 4 * 32M 32 110632783 110632814 This is working fine: counts = assays(summarizeOverlaps(tx, reads, mode = "Union"))$counts Method "IntersectionStrict" does not work: counts1 = assays(summarizeOverlaps(tx, liver, mode = "IntersectionStrict"))$counts Error in assays(summarizeOverlaps(tx, liver, mode = "IntersectionStrict")) : error in evaluating the argument 'x' in selecting a method for function 'assays': Error in queryHits(findOverlaps(query, subject, maxgap = maxgap, minoverlap = minoverlap, : error in evaluating the argument 'x' in selecting a method for function 'queryHits': Error in .findOverlaps.circle(circle.length, seqselect(queryRanges, qIdxs), : overlap type "within" is not yet supported for circular sequence MT What can I do? Do I have to omit all genes on the MT? Best, Stefanie

I forgot: > sessionInfo() R version 2.15.0 (2012-03-30) Platform: x86_64-unknown-linux-gnu (64-bit) locale: [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8 [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8 [7] LC_PAPER=C LC_NAME=C [9] LC_ADDRESS=C LC_TELEPHONE=C [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C attached base packages: [1] stats graphics grDevices utils datasets methods base other attached packages: [1] biomaRt_2.12.0 GenomicFeatures_1.8.1 AnnotationDbi_1.18.0 [4] Biobase_2.16.0 Rsamtools_1.8.0 Biostrings_2.24.0 [7] GenomicRanges_1.8.1 IRanges_1.14.2 BiocGenerics_0.2.0 loaded via a namespace (and not attached): [1] bitops_1.0-4.1 BSgenome_1.24.0 DBI_0.2-5 RCurl_1.91-1 [5] RSQLite_0.11.1 rtracklayer_1.16.0 stats4_2.15.0 tools_2.15.0 [9] XML_3.9-4 zlibbioc_1.2.0

Hi Stefanie, On 08/01/12 05:41, Stefanie Tauber wrote: > Dear list, > > when I try to count human reads vs a human transcript database, I get the > following error: > > > library(GenomicFeatures) > library(biomaRt) > > Here I construct my human transcript database: > > ensembl = useDataset("hsapiens_gene_ensembl",mart=useMart("ensembl")) > > res = getBM(attributes = "ensembl_transcript_id", filters = > ("biotype","status"),values=list("protein_coding","known"),mart = ensembl) > > humanDb = makeTranscriptDbFromBiomart(biomart = "ensembl", dataset = > "hsapiens_gene_ensembl",transcript_ids = as.character(res[,1])) > > tx = transcriptsBy(humanDb) > > > The reads are an GappedAlignments object: > reads > > GappedAlignments with 11791172 alignments and 0 elementMetadata cols: > seqnames strand cigar qwidth start end > <rle> <rle> <character> <integer> <integer> <integer> > SRR015293.3 3 * 32M 32 186338939 186338970 > SRR015293.5 16 * 32M 32 72094409 72094440 > SRR015293.7 1 * 32M 32 159683461 159683492 > SRR015293.9 4 * 32M 32 155529684 155529715 > SRR015293.10 4 * 32M 32 110632783 110632814 > > > This is working fine: > counts = assays(summarizeOverlaps(tx, reads, mode = "Union"))$counts > > Method "IntersectionStrict" does not work: > counts1 = assays(summarizeOverlaps(tx, liver, mode = "IntersectionStrict"))$counts > Error in assays(summarizeOverlaps(tx, liver, mode = "IntersectionStrict")) : > error in evaluating the argument 'x' in selecting a method for function > 'assays': Error in queryHits(findOverlaps(query, subject, maxgap = maxgap, > minoverlap = minoverlap, : > error in evaluating the argument 'x' in selecting a method for function > 'queryHits': Error in .findOverlaps.circle(circle.length, seqselect(queryRanges, > qIdxs), : > overlap type "within" is not yet supported for circular sequence MT > > > What can I do? Do I have to omit all genes on the MT? Thanks for reporting this. The problem was that findOverlaps for GRanges does not currently support type="within" for circular chromosomes. The mode 'IntersectionStrict' uses findOverlaps in this way whereas 'Union' and 'IntersectionNotEmpty' do not. I have added a check for circular chromosomes to the 'IntersectionStrict' mode. This will issue a warning when circular chromosomes are encountered and they will be omitted from the counting. This has been fixed in GenomicRanges release 1.8.11 and devel 1.9.41. These versions will be available through biocLite() tomorrow morning (~9 am PST) or through svn immediately. Valerie > > Best, > Stefanie > > _______________________________________________ > Bioconductor mailing list > Bioconductor at r-project.org > https://stat.ethz.ch/mailman/listinfo/bioconductor > Search the archives: http://news.gmane.org/gmane.science.biology.informatics.conductor