countBamInGRanges.exomeDepth | R Documentation |
Parses a BAM file and count reads that are located within a target region defined by a GenomicRanges object.
countBamInGRanges.exomeDepth( bam.file, index = bam.file, granges, min.mapq = 1, read.width = 1 )
bam.file |
BAM file to be parsed |
index |
Index of the BAM file, without the '.bai' suffix. |
granges |
Genomic ranges object defining the bins |
min.mapq |
Minimum read mapping quality (Phred scaled). |
read.width |
For single end reads, an estimate of the frament size. For paired reads, the fragment size can be directly computed from the paired alignment and this value is ignored. |
Largely derived from its equivalent function in the exomeCopy package.
A GRanges object with count data.
minimum_bam_file <- system.file('extdata/minimum_1_25630000_25650000.bam', package = 'ExomeDepth') data(exons.hg19) exons.hg19.RHD <- subset(exons.hg19, grepl(pattern = '^RHD', exons.hg19[['name']])) my_range <- GenomicRanges::GRanges(paste0(exons.hg19.RHD$chromosome, ":", exons.hg19.RHD$start, '-', exons.hg19.RHD$end)) print(countBamInGRanges.exomeDepth(bam.file = minimum_bam_file, granges = my_range))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.