This function takes bam file paths and read them into GRanges objects. Note: Can be quite lengthy for .bam files with 5+ millions fragments.

importPEBamFiles(
  files,
  genome = NULL,
  where = NULL,
  max_insert_size = 1000,
  shift_ATAC_fragments = FALSE,
  cores = 10,
  verbose = TRUE
)

Arguments

files

character vector, each element of the vector is the path of an individual .bam file.

genome

character, genome ID (e.g. sacCer3, ce11, dm6, danRer10, mm10 or hg38).

where

GRanges, only import the fragments mapping to the input GRanges (can fasten the import process a lot).

max_insert_size

Integer, filter out fragments larger than this size.

shift_ATAC_fragments

Boolean, if the fragments come from ATAC-seq, one might want to shift the extremities by +5 / -4 bp.

cores

Integer, number of cores to use when indexing bam files

verbose

Boolean

Value

A GRanges object containing fragments from the input .bam file.

Examples

bamfile <- system.file("extdata", "ex1.bam", package = "Rsamtools") fragments <- importPEBamFiles( bamfile, shift_ATAC_fragments = TRUE )
#> > Importing /Users/runner/work/_temp/Library/Rsamtools/extdata/ex1.bam ...
#> > Filtering /Users/runner/work/_temp/Library/Rsamtools/extdata/ex1.bam ...
#> > Shifting /Users/runner/work/_temp/Library/Rsamtools/extdata/ex1.bam ...
#> > /Users/runner/work/_temp/Library/Rsamtools/extdata/ex1.bam import completed.
fragments
#> GRanges object with 1572 ranges and 0 metadata columns: #> seqnames ranges strand #> <Rle> <IRanges> <Rle> #> [1] seq1 41-215 + #> [2] seq1 54-255 + #> [3] seq1 56-258 + #> [4] seq1 65-255 + #> [5] seq1 65-265 + #> ... ... ... ... #> [1568] seq2 1326-1542 - #> [1569] seq2 1336-1544 - #> [1570] seq2 1358-1550 - #> [1571] seq2 1380-1557 - #> [1572] seq2 1353-1562 - #> ------- #> seqinfo: 2 sequences from an unspecified genome; no seqlengths