Last data update: 2014.03.03

R: GAlignmentsList objects
GAlignmentsList-classR Documentation

GAlignmentsList objects

Description

The GAlignmentsList class is a container for storing a collection of GAlignments objects.

Details

A GAlignmentsList object contains a list of GAlignments objects. The majority of operations on this page are described in more detail on the GAlignments man page, see ?GAlignments.

Constructor

GAlignmentsList(...): Creates a GAlignmentsList from a list of GAlignments objects.

Accessors

In the code snippets below, x is a GAlignmentsList object.

length(x): Return the number of elements in x.

names(x), names(x) <- value: Get or set the names of the elements of x.

seqnames(x), seqnames(x) <- value: Get or set the name of the reference sequences of the alignments in each element of x.

rname(x), rname(x) <- value: Same as seqnames(x) and seqnames(x) <- value.

strand(x), strand(x) <- value: Get or set the strand of the alignments in each element of x.

cigar(x): Returns a character list of length length(x) containing the CIGAR string for the alignments in each element of x.

qwidth(x): Returns an integer list of length length(x) containing the length of the alignments in each element of x *after* hard clipping (i.e. the length of the query sequence that is stored in the corresponding SAM/BAM record).

start(x), end(x): Returns an integer list of length length(x) containing the "start" and "end" (respectively) of the alignments in each element of x.

width(x): Returns an integer list of length length(x) containing the "width" of the alignments in each element of x.

njunc(x): Returns an integer list of length x containing the number of junctions (i.e. N operations in the CIGAR) for the alignments in each element of x.

seqinfo(x), seqinfo(x) <- value: Get or set the information about the underlying sequences. value must be a Seqinfo object.

seqlevels(x), seqlevels(x) <- value: Get or set the sequence levels of the alignments in each element of x.

seqlengths(x), seqlengths(x) <- value: Get or set the sequence lengths for each element of x. seqlengths(x) is equivalent to seqlengths(seqinfo(x)). value can be a named non-negative integer or numeric vector eventually with NAs.

isCircular(x), isCircular(x) <- value: Get or set the circularity flags for the alignments in each element in x. value must be a named logical list eventually with NAs.

genome(x), genome(x) <- value: Get or set the genome identifier or assembly name for the alignments in each element of x. value must be a named character list eventually with NAs.

seqnameStyle(x): Get or set the seqname style for alignments in each element of x.

Coercion

In the code snippets below, x is a GAlignmentsList object.

granges(x, use.names=TRUE, use.mcols=FALSE, ignore.strand=FALSE), ranges(x, use.names=TRUE, use.mcols=FALSE):

Return either a GRanges or a IRanges object of length length(x). Note this coercion IGNORES the cigar information. The resulting ranges span the entire range, including any junctions or spaces between paired-end reads.

If use.names is TRUE, then the names on x (if any) are propagated to the returned object. If use.mcols is TRUE, then the metadata columns on x (if any) are propagated to the returned object.

granges coercion supports ignore.strand to allow ranges of opposite strand to be combined (see examples). All ranges in the resulting GRanges will have strand ‘*’.

grglist(x, use.names=TRUE, use.mcols=FALSE, ignore.strand=FALSE), rglist(x, use.names=TRUE, use.mcols=FALSE):

Return either a GRangesList or an IRangesList object of length length(x). This coercion RESPECTS the cigar information. The resulting ranges are fragments of the original ranges that do not include junctions or spaces between paired-end reads.

If use.names is TRUE, then the names on x (if any) are propagated to the returned object. If use.mcols is TRUE, then the metadata columns on x (if any) are propagated to the returned object.

grglist coercion supports ignore.strand to allow ranges of opposite strand to be combined (see examples). When ignore.strand is TRUE all ranges in the resulting GRangesList have strand ‘*’.

as(x, "GRanges"), as(x, "Ranges"), as(x, "GRangesList"), as(x, "RangesList"): Alternate ways of doing granges(x, use.names=TRUE, use.mcols=TRUE), ranges(x, use.names=TRUE, use.mcols=TRUE), grglist(x, use.names=TRUE, use.mcols=TRUE), and rglist(x, use.names=TRUE, use.mcols=TRUE), respectively.

as.data.frame(x, row.names = NULL, optional = FALSE, ..., value.name = "value", use.outer.mcols = FALSE, group_name.as.factor = FALSE): Coerces x to a data.frame. See as.data.frame on the List man page for details (?List).

as(x, "GALignmentsList"): Here x is a GAlignmentPairs object. Return a GAlignmentsList object of length length(x) where the i-th list element represents the ranges of the i-th alignment pair in x.

Subsetting and related operations

In the code snippets below, x is a GAlignmentsList object.

x[i], x[i] <- value: Get or set list elements i. i can be a numeric or logical vector. value must be a GAlignments.

x[[i]], x[[i]] <- value: Same as x[i], x[i] <- value.

x[i, j], x[i, j] <- value: Get or set list elements i with optional metadata columns j. i can be a numeric, logical or missing. value must be a GAlignments.

Combining

c(...): Concatenates the GAlignmentsList objects in ....

Author(s)

Valerie Obenchain <vobencha@fredhutch.org

References

http://samtools.sourceforge.net/

See Also

  • readGAlignmentsList for reading genomic alignments from a file (typically a BAM file) into a GAlignmentsList object.

  • GAlignments and GAlignmentPairs objects for handling aligned single- and paired-end reads, respectively.

  • junctions-methods for extracting and summarizing junctions from a GAlignmentsList object.

  • findOverlaps-methods for finding range overlaps between a GAlignmentsList object and another range-based object.

  • seqinfo in the GenomeInfoDb package for getting/setting/modifying the sequence information stored in an object.

  • The GRanges and GRangesList classes defined and documented in the GenomicRanges package.

Examples

gal1 <- GAlignments(
    seqnames=Rle(factor(c("chr1", "chr2", "chr1", "chr3")),
        c(1, 3, 2, 4)),
    pos=1:10, cigar=paste0(10:1, "M"),
    strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
    names=head(letters, 10), score=1:10)

gal2 <- GAlignments(
    seqnames=Rle(factor(c("chr2", "chr4")), c(3, 4)), pos=1:7,
    cigar=c("5M", "3M2N3M2N3M", "5M", "10M", "5M1N4M", "8M2N1M", "5M"),
    strand=Rle(strand(c("-", "+")), c(4, 3)),
    names=tail(letters, 7), score=1:7)

galist <- GAlignmentsList(noGaps=gal1, Gaps=gal2)


## ---------------------------------------------------------------------
## A. BASIC MANIPULATION
## ---------------------------------------------------------------------

length(galist)
names(galist)
seqnames(galist)
strand(galist)
head(cigar(galist))
head(qwidth(galist))
head(start(galist))
head(end(galist))
head(width(galist))
head(njunc(galist))
seqlevels(galist)

## Rename the reference sequences:
seqlevels(galist) <- sub("chr", "seq", seqlevels(galist))
seqlevels(galist)

grglist(galist)  # a GRangesList object
rglist(galist)   # an IRangesList object

## ---------------------------------------------------------------------
## B. SUBSETTING
## ---------------------------------------------------------------------

galist[strand(galist) == "-"]
has_junctions <- sapply(galist,
                        function(x) any(grepl("N", cigar(x), fixed=TRUE)))
galist[has_junctions]

## Different ways to subset:
galist[2]             # a GAlignments object of length 1
galist[[2]]           # a GAlignments object of length 1
grglist(galist[2])  # a GRangesList object of length 1
rglist(galist[2])   # a NormalIRangesList object of length 1

## ---------------------------------------------------------------------
## C. mcols()/elementMetadata()
## ---------------------------------------------------------------------

## Metadata can be defined on the individual GAlignment elements
## and the overall GAlignmentsList object. By default, 'level=between' 
## extracts the GALignmentsList metadata. Using 'level=within' 
## will extract the metadata on the individual GAlignments objects.

mcols(galist) ## no metadata on the GAlignmentsList object
mcols(galist, level="within")


## ---------------------------------------------------------------------
## D. readGAlignmentsList()
## ---------------------------------------------------------------------

library(pasillaBamSubset)

## 'file' as character.
fl <- untreated3_chr4()
galist1 <- readGAlignmentsList(fl)

galist1[1:3]
length(galist1)
table(elementNROWS(galist1))

## When 'file' is a BamFile, 'asMates' must be TRUE. If FALSE,
## the data are treated as single-end and each list element of the
## GAlignmentsList will be of length 1. For single-end data 
## use readGAlignments() instead of readGAlignmentsList().
bf <- BamFile(fl, yieldSize=3, asMates=TRUE)
readGAlignmentsList(bf)

## Use a 'param' to fine tune the results.
param <- ScanBamParam(flag=scanBamFlag(isProperPair=TRUE))
galist2 <- readGAlignmentsList(fl, param=param)
length(galist2)


## ---------------------------------------------------------------------
## E. COERCION 
## ---------------------------------------------------------------------

## The granges() and grlist() coercions support 'ignore.strand' to 
## allow ranges from different strands to be combined. In this example 
## paired-end reads aligned to opposite strands were read into a 
## GAlignmentsList. If the desired operation is to combine these ranges, 
## regardless of junctions or the space between pairs, 'ignore.strand'
## must be TRUE.
granges(galist[1])
granges(galist[1], ignore.strand=TRUE)

## grglist()
galist <- GAlignmentsList(noGaps=gal1, Gaps=gal2)
grglist(galist)
grglist(galist, ignore.strand=TRUE)

Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(GenomicAlignments)
Loading required package: BiocGenerics
Loading required package: parallel

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:parallel':

    clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
    clusterExport, clusterMap, parApply, parCapply, parLapply,
    parLapplyLB, parRapply, parSapply, parSapplyLB

The following objects are masked from 'package:stats':

    IQR, mad, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, append,
    as.data.frame, cbind, colnames, do.call, duplicated, eval, evalq,
    get, grep, grepl, intersect, is.unsorted, lapply, lengths, mapply,
    match, mget, order, paste, pmax, pmax.int, pmin, pmin.int, rank,
    rbind, rownames, sapply, setdiff, sort, table, tapply, union,
    unique, unsplit

Loading required package: S4Vectors
Loading required package: stats4

Attaching package: 'S4Vectors'

The following objects are masked from 'package:base':

    colMeans, colSums, expand.grid, rowMeans, rowSums

Loading required package: IRanges
Loading required package: GenomeInfoDb
Loading required package: GenomicRanges
Loading required package: SummarizedExperiment
Loading required package: Biobase
Welcome to Bioconductor

    Vignettes contain introductory material; view with
    'browseVignettes()'. To cite Bioconductor, see
    'citation("Biobase")', and for packages 'citation("pkgname")'.

Loading required package: Biostrings
Loading required package: XVector
Loading required package: Rsamtools
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/GenomicAlignments/GAlignmentsList-class.Rd_%03d_medium.png", width=480, height=480)
> ### Name: GAlignmentsList-class
> ### Title: GAlignmentsList objects
> ### Aliases: class:GAlignmentsList GAlignmentsList-class GAlignmentsList
> ###   GAlignmentsList updateObject,GAlignmentsList-method
> ###   names,GAlignmentsList-method names<-,GAlignmentsList-method
> ###   seqnames,GAlignmentsList-method seqnames<-,GAlignmentsList-method
> ###   rname,GAlignmentsList-method rname<-,GAlignmentsList-method
> ###   strand,GAlignmentsList-method strand<-,GAlignmentsList-method
> ###   strand<-,GAlignmentsList,character-method
> ###   cigar,GAlignmentsList-method qwidth,GAlignmentsList-method
> ###   njunc,GAlignmentsList-method elementMetadata,GAlignmentsList-method
> ###   elementMetadata<-,GAlignmentsList-method
> ###   seqinfo,GAlignmentsList-method seqinfo<-,GAlignmentsList-method
> ###   start,GAlignmentsList-method end,GAlignmentsList-method
> ###   width,GAlignmentsList-method ranges,GAlignmentsList-method
> ###   granges,GAlignmentsList-method grglist,GAlignmentsList-method
> ###   rglist,GAlignmentsList-method coerce,GAlignmentsList,Ranges-method
> ###   coerce,GAlignmentsList,GRanges-method
> ###   coerce,GAlignmentsList,GRangesList-method
> ###   coerce,GAlignmentsList,RangesList-method
> ###   coerce,GAlignmentPairs,GAlignmentsList-method
> ###   c,GAlignmentsList-method relistToClass,GAlignments-method
> ###   show,GAlignmentsList-method
> ### Keywords: methods classes
> 
> ### ** Examples
> 
> gal1 <- GAlignments(
+     seqnames=Rle(factor(c("chr1", "chr2", "chr1", "chr3")),
+         c(1, 3, 2, 4)),
+     pos=1:10, cigar=paste0(10:1, "M"),
+     strand=Rle(strand(c("-", "+", "*", "+", "-")), c(1, 2, 2, 3, 2)),
+     names=head(letters, 10), score=1:10)
> 
> gal2 <- GAlignments(
+     seqnames=Rle(factor(c("chr2", "chr4")), c(3, 4)), pos=1:7,
+     cigar=c("5M", "3M2N3M2N3M", "5M", "10M", "5M1N4M", "8M2N1M", "5M"),
+     strand=Rle(strand(c("-", "+")), c(4, 3)),
+     names=tail(letters, 7), score=1:7)
> 
> galist <- GAlignmentsList(noGaps=gal1, Gaps=gal2)
> 
> 
> ## ---------------------------------------------------------------------
> ## A. BASIC MANIPULATION
> ## ---------------------------------------------------------------------
> 
> length(galist)
[1] 2
> names(galist)
[1] "noGaps" "Gaps"  
> seqnames(galist)
RleList of length 2
$noGaps
factor-Rle of length 10 with 4 runs
  Lengths:    1    3    2    4
  Values : chr1 chr2 chr1 chr3
Levels(4): chr1 chr2 chr3 chr4

$Gaps
factor-Rle of length 7 with 2 runs
  Lengths:    3    4
  Values : chr2 chr4
Levels(4): chr1 chr2 chr3 chr4

> strand(galist)
RleList of length 2
$noGaps
factor-Rle of length 10 with 5 runs
  Lengths: 1 2 2 3 2
  Values : - + * + -
Levels(3): + - *

$Gaps
factor-Rle of length 7 with 2 runs
  Lengths: 4 3
  Values : - +
Levels(3): + - *

> head(cigar(galist))
CharacterList of length 2
[["noGaps"]] 10M 9M 8M 7M 6M 5M 4M 3M 2M 1M
[["Gaps"]] 5M 3M2N3M2N3M 5M 10M 5M1N4M 8M2N1M 5M
> head(qwidth(galist))
IntegerList of length 2
[["noGaps"]] 10 9 8 7 6 5 4 3 2 1
[["Gaps"]] 5 9 5 10 9 9 5
> head(start(galist))
IntegerList of length 2
[["noGaps"]] 1 2 3 4 5 6 7 8 9 10
[["Gaps"]] 1 2 3 4 5 6 7
> head(end(galist))
IntegerList of length 2
[["noGaps"]] 10 10 10 10 10 10 10 10 10 10
[["Gaps"]] 5 14 7 13 14 16 11
> head(width(galist))
IntegerList of length 2
[["noGaps"]] 10 9 8 7 6 5 4 3 2 1
[["Gaps"]] 5 13 5 10 10 11 5
> head(njunc(galist))
IntegerList of length 2
[["noGaps"]] 0 0 0 0 0 0 0 0 0 0
[["Gaps"]] 0 2 0 0 1 1 0
> seqlevels(galist)
[1] "chr1" "chr2" "chr3" "chr4"
> 
> ## Rename the reference sequences:
> seqlevels(galist) <- sub("chr", "seq", seqlevels(galist))
> seqlevels(galist)
[1] "seq1" "seq2" "seq3" "seq4"
> 
> grglist(galist)  # a GRangesList object
GRangesList object of length 2:
$noGaps 
GRanges object with 10 ranges and 0 metadata columns:
       seqnames    ranges strand
          <Rle> <IRanges>  <Rle>
   [1]     seq1  [ 1, 10]      -
   [2]     seq2  [ 2, 10]      +
   [3]     seq2  [ 3, 10]      +
   [4]     seq2  [ 4, 10]      *
   [5]     seq1  [ 5, 10]      *
   [6]     seq1  [ 6, 10]      +
   [7]     seq3  [ 7, 10]      +
   [8]     seq3  [ 8, 10]      +
   [9]     seq3  [ 9, 10]      -
  [10]     seq3  [10, 10]      -

...
<1 more element>
-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> rglist(galist)   # an IRangesList object
IRangesList of length 2
$noGaps
IRanges object with 10 ranges and 0 metadata columns:
           start       end     width
       <integer> <integer> <integer>
   [1]         1        10        10
   [2]         2        10         9
   [3]         3        10         8
   [4]         4        10         7
   [5]         5        10         6
   [6]         6        10         5
   [7]         7        10         4
   [8]         8        10         3
   [9]         9        10         2
  [10]        10        10         1

...
<1 more element>
> 
> ## ---------------------------------------------------------------------
> ## B. SUBSETTING
> ## ---------------------------------------------------------------------
> 
> galist[strand(galist) == "-"]
GAlignmentsList object of length 2:
$noGaps 
GAlignments object with 3 alignments and 1 metadata column:
    seqnames strand cigar qwidth start end width njunc | score
  a     seq1      -   10M     10     1  10    10     0 |     1
  i     seq3      -    2M      2     9  10     2     0 |     9
  j     seq3      -    1M      1    10  10     1     0 |    10

$Gaps 
GAlignments object with 4 alignments and 1 metadata column:
    seqnames strand      cigar qwidth start end width njunc | score
  t     seq2      -         5M      5     1   5     5     0 |     1
  u     seq2      - 3M2N3M2N3M      9     2  14    13     2 |     2
  v     seq2      -         5M      5     3   7     5     0 |     3
  w     seq4      -        10M     10     4  13    10     0 |     4

-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> has_junctions <- sapply(galist,
+                         function(x) any(grepl("N", cigar(x), fixed=TRUE)))
> galist[has_junctions]
GAlignmentsList object of length 1:
$Gaps 
GAlignments object with 7 alignments and 1 metadata column:
    seqnames strand      cigar qwidth start end width njunc | score
  t     seq2      -         5M      5     1   5     5     0 |     1
  u     seq2      - 3M2N3M2N3M      9     2  14    13     2 |     2
  v     seq2      -         5M      5     3   7     5     0 |     3
  w     seq4      -        10M     10     4  13    10     0 |     4
  x     seq4      +     5M1N4M      9     5  14    10     1 |     5
  y     seq4      +     8M2N1M      9     6  16    11     1 |     6
  z     seq4      +         5M      5     7  11     5     0 |     7

-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> 
> ## Different ways to subset:
> galist[2]             # a GAlignments object of length 1
GAlignmentsList object of length 1:
$Gaps 
GAlignments object with 7 alignments and 1 metadata column:
    seqnames strand      cigar qwidth start end width njunc | score
  t     seq2      -         5M      5     1   5     5     0 |     1
  u     seq2      - 3M2N3M2N3M      9     2  14    13     2 |     2
  v     seq2      -         5M      5     3   7     5     0 |     3
  w     seq4      -        10M     10     4  13    10     0 |     4
  x     seq4      +     5M1N4M      9     5  14    10     1 |     5
  y     seq4      +     8M2N1M      9     6  16    11     1 |     6
  z     seq4      +         5M      5     7  11     5     0 |     7

-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> galist[[2]]           # a GAlignments object of length 1
GAlignments object with 7 alignments and 1 metadata column:
    seqnames strand       cigar    qwidth     start       end     width
       <Rle>  <Rle> <character> <integer> <integer> <integer> <integer>
  t     seq2      -          5M         5         1         5         5
  u     seq2      -  3M2N3M2N3M         9         2        14        13
  v     seq2      -          5M         5         3         7         5
  w     seq4      -         10M        10         4        13        10
  x     seq4      +      5M1N4M         9         5        14        10
  y     seq4      +      8M2N1M         9         6        16        11
  z     seq4      +          5M         5         7        11         5
        njunc |     score
    <integer> | <integer>
  t         0 |         1
  u         2 |         2
  v         0 |         3
  w         0 |         4
  x         1 |         5
  y         1 |         6
  z         0 |         7
  -------
  seqinfo: 4 sequences from an unspecified genome; no seqlengths
> grglist(galist[2])  # a GRangesList object of length 1
GRangesList object of length 1:
$Gaps 
GRanges object with 11 ranges and 0 metadata columns:
       seqnames    ranges strand
          <Rle> <IRanges>  <Rle>
   [1]     seq2  [ 1,  5]      -
   [2]     seq2  [ 2,  4]      -
   [3]     seq2  [ 7,  9]      -
   [4]     seq2  [12, 14]      -
   [5]     seq2  [ 3,  7]      -
   [6]     seq4  [ 4, 13]      -
   [7]     seq4  [ 5,  9]      +
   [8]     seq4  [11, 14]      +
   [9]     seq4  [ 6, 13]      +
  [10]     seq4  [16, 16]      +
  [11]     seq4  [ 7, 11]      +

-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> rglist(galist[2])   # a NormalIRangesList object of length 1
IRangesList of length 1
$Gaps
IRanges object with 11 ranges and 0 metadata columns:
           start       end     width
       <integer> <integer> <integer>
   [1]         1         5         5
   [2]         2         4         3
   [3]         7         9         3
   [4]        12        14         3
   [5]         3         7         5
   [6]         4        13        10
   [7]         5         9         5
   [8]        11        14         4
   [9]         6        13         8
  [10]        16        16         1
  [11]         7        11         5

> 
> ## ---------------------------------------------------------------------
> ## C. mcols()/elementMetadata()
> ## ---------------------------------------------------------------------
> 
> ## Metadata can be defined on the individual GAlignment elements
> ## and the overall GAlignmentsList object. By default, 'level=between' 
> ## extracts the GALignmentsList metadata. Using 'level=within' 
> ## will extract the metadata on the individual GAlignments objects.
> 
> mcols(galist) ## no metadata on the GAlignmentsList object
DataFrame with 2 rows and 0 columns
> mcols(galist, level="within")
SplitDataFrameList of length 2
$noGaps
DataFrame with 10 rows and 1 column
       score
   <integer>
1          1
2          2
3          3
4          4
5          5
6          6
7          7
8          8
9          9
10        10

$Gaps
DataFrame with 7 rows and 1 column
      score
  <integer>
1         1
2         2
3         3
4         4
5         5
6         6
7         7

> 
> 
> ## ---------------------------------------------------------------------
> ## D. readGAlignmentsList()
> ## ---------------------------------------------------------------------
> 
> library(pasillaBamSubset)
> 
> ## 'file' as character.
> fl <- untreated3_chr4()
> galist1 <- readGAlignmentsList(fl)
> 
> galist1[1:3]
GAlignmentsList object of length 3:
[[1]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start end width njunc
  [1]     chr4      +   37M     37   169 205    37     0
  [2]     chr4      -   37M     37   326 362    37     0

[[2]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start  end width njunc
  [1]     chr4      +   37M     37   946  982    37     0
  [2]     chr4      -   37M     37   986 1022    37     0

[[3]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start  end width njunc
  [1]     chr4      +   37M     37   943  979    37     0
  [2]     chr4      -   37M     37  1086 1122    37     0

-------
seqinfo: 8 sequences from an unspecified genome
> length(galist1)
[1] 96636
> table(elementNROWS(galist1))

    1     2     3     4     5     6     7     8     9 
18191 78297    69    60     7     8     2     1     1 
> 
> ## When 'file' is a BamFile, 'asMates' must be TRUE. If FALSE,
> ## the data are treated as single-end and each list element of the
> ## GAlignmentsList will be of length 1. For single-end data 
> ## use readGAlignments() instead of readGAlignmentsList().
> bf <- BamFile(fl, yieldSize=3, asMates=TRUE)
> readGAlignmentsList(bf)
GAlignmentsList object of length 3:
[[1]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start end width njunc
  [1]     chr4      +   37M     37   169 205    37     0
  [2]     chr4      -   37M     37   326 362    37     0

[[2]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start  end width njunc
  [1]     chr4      +   37M     37   946  982    37     0
  [2]     chr4      -   37M     37   986 1022    37     0

[[3]] 
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand cigar qwidth start  end width njunc
  [1]     chr4      +   37M     37   943  979    37     0
  [2]     chr4      -   37M     37  1086 1122    37     0

-------
seqinfo: 8 sequences from an unspecified genome
> 
> ## Use a 'param' to fine tune the results.
> param <- ScanBamParam(flag=scanBamFlag(isProperPair=TRUE))
> galist2 <- readGAlignmentsList(fl, param=param)
> length(galist2)
[1] 45828
> 
> 
> ## ---------------------------------------------------------------------
> ## E. COERCION 
> ## ---------------------------------------------------------------------
> 
> ## The granges() and grlist() coercions support 'ignore.strand' to 
> ## allow ranges from different strands to be combined. In this example 
> ## paired-end reads aligned to opposite strands were read into a 
> ## GAlignmentsList. If the desired operation is to combine these ranges, 
> ## regardless of junctions or the space between pairs, 'ignore.strand'
> ## must be TRUE.
> granges(galist[1])
GRanges object with 7 ranges and 0 metadata columns:
         seqnames    ranges strand
            <Rle> <IRanges>  <Rle>
  noGaps     seq1   [6, 10]      +
  noGaps     seq1   [1, 10]      -
  noGaps     seq1   [5, 10]      *
  noGaps     seq2   [2, 10]      +
  noGaps     seq2   [4, 10]      *
  noGaps     seq3   [7, 10]      +
  noGaps     seq3   [9, 10]      -
  -------
  seqinfo: 4 sequences from an unspecified genome; no seqlengths
Warning message:
In .local(x, use.names, use.mcols, ...) :
  For some list elements in 'x', the ranges are not on the same chromosome and strand. Cannot extract a single range for them. As a consequence, the returned GRanges object is not parallel to 'x'. Consider using 'ignore.strand=TRUE'.
> granges(galist[1], ignore.strand=TRUE)
GRanges object with 3 ranges and 0 metadata columns:
         seqnames    ranges strand
            <Rle> <IRanges>  <Rle>
  noGaps     seq1   [1, 10]      *
  noGaps     seq2   [2, 10]      *
  noGaps     seq3   [7, 10]      *
  -------
  seqinfo: 4 sequences from an unspecified genome; no seqlengths
Warning message:
In .local(x, use.names, use.mcols, ...) :
  For some list elements in 'x', the ranges are not on the same chromosome and strand. Cannot extract a single range for them. As a consequence, the returned GRanges object is not parallel to 'x'.
> 
> ## grglist()
> galist <- GAlignmentsList(noGaps=gal1, Gaps=gal2)
> grglist(galist)
GRangesList object of length 2:
$noGaps 
GRanges object with 10 ranges and 0 metadata columns:
       seqnames    ranges strand
          <Rle> <IRanges>  <Rle>
   [1]     chr1  [ 1, 10]      -
   [2]     chr2  [ 2, 10]      +
   [3]     chr2  [ 3, 10]      +
   [4]     chr2  [ 4, 10]      *
   [5]     chr1  [ 5, 10]      *
   [6]     chr1  [ 6, 10]      +
   [7]     chr3  [ 7, 10]      +
   [8]     chr3  [ 8, 10]      +
   [9]     chr3  [ 9, 10]      -
  [10]     chr3  [10, 10]      -

...
<1 more element>
-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> grglist(galist, ignore.strand=TRUE)
GRangesList object of length 2:
$noGaps 
GRanges object with 10 ranges and 0 metadata columns:
       seqnames    ranges strand
          <Rle> <IRanges>  <Rle>
   [1]     chr1  [ 1, 10]      *
   [2]     chr2  [ 2, 10]      *
   [3]     chr2  [ 3, 10]      *
   [4]     chr2  [ 4, 10]      *
   [5]     chr1  [ 5, 10]      *
   [6]     chr1  [ 6, 10]      *
   [7]     chr3  [ 7, 10]      *
   [8]     chr3  [ 8, 10]      *
   [9]     chr3  [ 9, 10]      *
  [10]     chr3  [10, 10]      *

...
<1 more element>
-------
seqinfo: 4 sequences from an unspecified genome; no seqlengths
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>