Get fragment widths

Share:

Description

From a BAM file and a particular transcript (recommened to be the single isoform of a gene), this function returns estimates of the fragment widths, by mapping the fragment alignments to the transcript coordinates.

Usage

1
getFragmentWidths(bam.file, tx)

Arguments

bam.file

a character string pointing to a BAM file

tx

a GRanges object of the exons of a single isoform gene

Value

a numeric vector of estimated fragment widths

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
# these next lines just write out a BAM file from R
# typically you would already have a BAM file
library(alpineData)
library(GenomicAlignments)
library(rtracklayer)
gap <- ERR188088()
dir <- system.file(package="alpineData", "extdata")
bam.file <- c("ERR188088" = file.path(dir,"ERR188088.bam"))
export(gap, con=bam.file)

data(preprocessedData)

w <- getFragmentWidths(bam.file, ebt.fit[[2]])
quantile(w, c(.025, .975))