GAlignmentPairs objects

Description

The GAlignmentPairs class is a container for genomic alignment pairs.

IMPORTANT NOTE: The GAlignmentPairs container only supports pairs where the 2 alignments are on opposite strands of the same chromosome at the moment.

Details

A GAlignmentPairs object is a list-like object where each element describes a pair of genomic alignment.

An "alignment pair" is made of a "first" and a "last"/"second" alignment, and is formally represented by a GAlignments object of length 2. It is typically representing a hit of a paired-end read to the reference genome that was used by the aligner. More precisely, in a given pair, the "first" alignment represents the hit of the first end of the read (aka "first segment in the template", using SAM Spec terminology), and the "last" alignment represents the hit of the second end of the read (aka "last segment in the template", using SAM Spec terminology).

In general, a GAlignmentPairs object will be created by loading records from a BAM (or SAM) file containing aligned paired-end reads, using the readGAlignmentPairs function (see below). Each element in the returned object will be obtained by pairing 2 records.

Constructor

GAlignmentPairs(first, last, strandMode=1, isProperPair=TRUE, names=NULL): Low-level GAlignmentPairs constructor. Generally not used directly.

Accessors

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

strandMode(x), strandMode(x) <- value: The strand mode is a per-object switch on GAlignmentPairs objects that controls the behavior of the strand getter. More precisely, it indicates how the strand of a pair should be inferred from the strand of the first and last alignments in the pair:

• 0: strand of the pair is always *.

• 1: strand of the pair is strand of its first alignment. This mode should be used when the paired-end data was generated using one of the following stranded protocols: Directional Illumina (Ligation), Standard SOLiD.

• 2: strand of the pair is strand of its last alignment. This mode should be used when the paired-end data was generated using one of the following stranded protocols: dUTP, NSR, NNSR, Illumina stranded TruSeq PE protocol.

These modes are equivalent to strandSpecific equal 0, 1, and 2, respectively, for the featureCounts function defined in the Rsubread package.

Note that, by default, the readGAlignmentPairs function sets the strand mode to 1 on the returned GAlignmentPairs object. The function has a strandMode argument to let the user set a different strand mode. The strand mode can also be changed any time with the strandMode setter.

Also note that 3rd party programs TopHat2 and Cufflinks have a --library-type option to let the user specify which protocol was used. Please refer to the documentation of these programs for more information.

length(x): Return the number of alignment pairs in x.

names(x), names(x) <- value: Get or set the names on x. See readGAlignmentPairs for how to automatically extract and set the names when reading the alignments from a file.

first(x, real.strand=FALSE), last(x, real.strand=FALSE): second(x, real.strand=FALSE): Get the "first" or "last"/"second" alignment for each alignment pair in x. The result is a GAlignments object of the same length as x.

If real.strand=TRUE, then the strand is inverted on-the-fly according to the strand mode currently set on the object (see strandMode(x) above). More precisely, if strandMode(x) is 0, then the strand is set to * for the GAlignments object returned by both, first() and last(). If strandMode(x) is 1, then the strand of the object returned by last() is inverted. If strandMode(x) is 2, then the strand of the object returned by first() is inverted.

seqnames(x): Get the name of the reference sequence for each alignment pair in x. When reading the alignments from a BAM file, this comes from the RNAME field which has the same value for the 2 records in a pair (readGAlignmentPairs, the function used for reading paired-end reads from a BAM file as a GAlignmentPairs object, rejects pairs with incompatible RNAME values).

strand(x), strand(x) <- value: Get or set the strand for each alignment pair in x. Obeys strandMode(x) above to infer the strand of a pair.

njunc(x): Equivalent to njunc(first(x)) + njunc(last(x)).

isProperPair(x): Get the "isProperPair" flag bit (bit 0x2 in SAM Spec) set by the aligner for each alignment pair in 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. seqlevels(x) is equivalent to seqlevels(seqinfo(x)) or to levels(seqnames(x)), those 2 expressions being guaranteed to return identical character vectors on a GAlignmentPairs object. value must be a character vector with no NAs. See ?seqlevels for more information.

seqlengths(x), seqlengths(x) <- value: Get or set the sequence lengths. 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. isCircular(x) is equivalent to isCircular(seqinfo(x)). value must be a named logical vector eventually with NAs.

genome(x), genome(x) <- value: Get or set the genome identifier or assembly name for each sequence. genome(x) is equivalent to genome(seqinfo(x)). value must be a named character vector eventually with NAs.

seqnameStyle(x): Get or set the seqname style for x. Note that this information is not stored in x but inferred by looking up seqnames(x) against a seqname style database stored in the seqnames.db metadata package (required). seqnameStyle(x) is equivalent to seqnameStyle(seqinfo(x)) and can return more than 1 seqname style (with a warning) in case the style cannot be determined unambiguously.

Vector methods

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

x[i]: Return a new GAlignmentPairs object made of the selected alignment pairs.

List methods

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

x[[i]]: Extract the i-th alignment pair as a GAlignments object of length 2. As expected x[[i]][1] and x[[i]][2] are respectively the "first" and "last" alignments in the pair.

unlist(x, use.names=TRUE): Return the GAlignments object conceptually defined by c(x[[1]], x[[2]], ..., x[[length(x)]]). use.names determines whether x names should be propagated to the result or not.

Coercion

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

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

Return a GRanges object (for granges()) or IRanges) object (for ranges()) parallel to x where the i-th element is the range of the genomic region spanned by the i-th alignment in x. All gaps in the region are ignored.

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(x, use.mcols=FALSE, drop.D.ranges=FALSE):

Return a GRangesList object of length length(x) where the i-th element represents the ranges (with respect to the reference) of the i-th alignment pair in x. The strand of the returned ranges obeys the strand mode currently set on the object (see strandMode(x) above).

More precisely, if grl1 and grl2 are grglist(first(x, real.strand=TRUE), order.as.in.query=TRUE) and grglist(last(x, real.strand=TRUE), order.as.in.query=TRUE), respectively, then the i-th element in the returned GRangesList object is c(grl1[[i]], grl2[[i]]), if strandMode(x) is 1, or c(grl2[[i]], grl1[[i]]), if strandMode(x) is 2.

Note that this results in the ranges being always ordered consistently with the original "query template", that is, being in the order defined by walking the "query template" from the beginning to the end.

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.

If drop.D.ranges is TRUE, then deletions (Ds in the CIGAR) are treated like junctions (Ns in the CIGAR), that is, the ranges corresponding to deletions are dropped.

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

as(x, "GAlignments"): Equivalent of unlist(x, use.names=TRUE).

Other methods

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

show(x): By default the show method displays 5 head and 5 tail elements. This can be changed by setting the global options showHeadLines and showTailLines. If the object length is less than (or equal to) the sum of these 2 options plus 1, then the full object is displayed. Note that these options also affect the display of GRanges and GAlignments objects, as well as other objects defined in the IRanges and Biostrings packages (e.g. Ranges and XStringSet objects).

Author(s)

Herv<c3><83><c2><a9> Pag<c3><83><c2><a8>s

See Also

• readGAlignmentPairs for reading aligned paired-end reads from a file (typically a BAM file) into a GAlignmentPairs object.

• GAlignments objects for handling aligned single-end reads.

• makeGAlignmentPairs for pairing the elements of a GAlignments object into a GAlignmentPairs object.

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

• coverage-methods for computing the coverage of a GAlignmentPairs object.

• findOverlaps-methods for finding range overlaps between a GAlignmentPairs 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

library(Rsamtools)  # for the ex1.bam file
ex1_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
galp <- readGAlignmentPairs(ex1_file, use.names=TRUE, strandMode=1)
galp

length(galp)
head(galp)
head(names(galp))

first(galp)
last(galp)
# or
second(galp)

strandMode(galp)
first(galp, real.strand=TRUE)
last(galp, real.strand=TRUE)
strand(galp)

strandMode(galp) <- 2
first(galp, real.strand=TRUE)
last(galp, real.strand=TRUE)
strand(galp)

seqnames(galp)

head(njunc(galp))
table(isProperPair(galp))
seqlevels(galp)

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

galp[[1]]
unlist(galp)

grglist(galp)  # a GRangesList object

strandMode(galp) <- 1
grglist(galp)

stopifnot(identical(unname(elementNROWS(grglist(galp))), njunc(galp) + 2L))

granges(galp)  # a GRanges object

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/GAlignmentPairs-class.Rd_%03d_medium.png", width=480, height=480)
> ### Name: GAlignmentPairs-class
> ### Title: GAlignmentPairs objects
> ### Aliases: class:GAlignmentPairs GAlignmentPairs-class GAlignmentPairs
> ###   strandMode strandMode,GAlignmentPairs-method strandMode<-
> ###   strandMode<-,GAlignmentPairs-method length,GAlignmentPairs-method
> ###   names,GAlignmentPairs-method names<-,GAlignmentPairs-method first
> ###   first,GAlignmentPairs-method second second,GAlignmentPairs-method
> ###   last last,GAlignmentPairs-method seqnames,GAlignmentPairs-method
> ###   strand,GAlignmentPairs-method strand<-,GAlignmentPairs-method
> ###   njunc,GAlignmentPairs-method isProperPair
> ###   isProperPair,GAlignmentPairs-method seqinfo,GAlignmentPairs-method
> ###   seqlevelsInUse,GAlignmentPairs-method
> ###   seqinfo<-,GAlignmentPairs-method [[,GAlignmentPairs,ANY,ANY-method
> ###   unlist,GAlignmentPairs-method ranges,GAlignmentPairs-method
> ###   granges,GAlignmentPairs-method grglist,GAlignmentPairs-method
> ###   coerce,GAlignmentPairs,Ranges-method
> ###   coerce,GAlignmentPairs,GRanges-method
> ###   coerce,GAlignmentPairs,GRangesList-method
> ###   coerce,GAlignmentPairs,GAlignments-method show,GAlignmentPairs-method
> ###   c,GAlignmentPairs-method left left,GAlignmentPairs-method right
> ###   right,GAlignmentPairs-method
> ### Keywords: methods classes
> 
> ### ** Examples
> 
> library(Rsamtools)  # for the ex1.bam file
> ex1_file <- system.file("extdata", "ex1.bam", package="Rsamtools")
> galp <- readGAlignmentPairs(ex1_file, use.names=TRUE, strandMode=1)
> galp
GAlignmentPairs object with 1572 pairs, strandMode=1, and 0 metadata columns:
                           seqnames strand   :       ranges  --       ranges
                              <Rle>  <Rle>   :    <IRanges>  --    <IRanges>
     EAS54_61:4:143:69:578     seq1      +   :     [36, 70]  --   [185, 219]
      B7_593:4:106:316:452     seq1      +   :     [49, 84]  --   [224, 259]
    EAS54_65:3:321:311:983     seq1      +   :     [51, 85]  --   [228, 262]
       B7_591:5:42:540:501     seq1      +   :     [60, 95]  --   [224, 259]
    EAS192_3:5:223:142:410     seq1      +   :     [60, 94]  --   [235, 269]
                       ...      ...    ... ...          ... ...          ...
  EAS139_11:5:52:1278:1478     seq2      -   : [1513, 1547]  -- [1322, 1356]
     EAS1_97:4:274:287:423     seq2      -   : [1515, 1549]  -- [1332, 1366]
    EAS54_71:8:105:854:975     seq2      -   : [1523, 1555]  -- [1354, 1388]
  EAS139_11:7:50:1229:1313     seq2      -   : [1528, 1562]  -- [1376, 1410]
     EAS114_26:7:37:79:581     seq2      -   : [1533, 1567]  -- [1349, 1383]
  -------
  seqinfo: 2 sequences from an unspecified genome
> 
> length(galp)
[1] 1572
> head(galp)
GAlignmentPairs object with 6 pairs, strandMode=1, and 0 metadata columns:
                         seqnames strand :    ranges --     ranges
                            <Rle>  <Rle> : <IRanges> --  <IRanges>
   EAS54_61:4:143:69:578     seq1      + :  [36, 70] -- [185, 219]
    B7_593:4:106:316:452     seq1      + :  [49, 84] -- [224, 259]
  EAS54_65:3:321:311:983     seq1      + :  [51, 85] -- [228, 262]
     B7_591:5:42:540:501     seq1      + :  [60, 95] -- [224, 259]
  EAS192_3:5:223:142:410     seq1      + :  [60, 94] -- [235, 269]
  EAS56_61:6:227:259:597     seq1      + :  [61, 95] -- [248, 282]
  -------
  seqinfo: 2 sequences from an unspecified genome
> head(names(galp))
[1] "EAS54_61:4:143:69:578"  "B7_593:4:106:316:452"   "EAS54_65:3:321:311:983"
[4] "B7_591:5:42:540:501"    "EAS192_3:5:223:142:410" "EAS56_61:6:227:259:597"
> 
> first(galp)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      +         35M        35        36
      B7_593:4:106:316:452     seq1      +         36M        36        49
    EAS54_65:3:321:311:983     seq1      +         35M        35        51
       B7_591:5:42:540:501     seq1      +         36M        36        60
    EAS192_3:5:223:142:410     seq1      +         35M        35        60
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      -         35M        35      1513
     EAS1_97:4:274:287:423     seq2      -         35M        35      1515
    EAS54_71:8:105:854:975     seq2      -         33M        33      1523
  EAS139_11:7:50:1229:1313     seq2      -         35M        35      1528
     EAS114_26:7:37:79:581     seq2      -         35M        35      1533
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578        70        35         0
      B7_593:4:106:316:452        84        36         0
    EAS54_65:3:321:311:983        85        35         0
       B7_591:5:42:540:501        95        36         0
    EAS192_3:5:223:142:410        94        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1547        35         0
     EAS1_97:4:274:287:423      1549        35         0
    EAS54_71:8:105:854:975      1555        33         0
  EAS139_11:7:50:1229:1313      1562        35         0
     EAS114_26:7:37:79:581      1567        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> last(galp)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      -         35M        35       185
      B7_593:4:106:316:452     seq1      -         36M        36       224
    EAS54_65:3:321:311:983     seq1      -         35M        35       228
       B7_591:5:42:540:501     seq1      -         36M        36       224
    EAS192_3:5:223:142:410     seq1      -         35M        35       235
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      +         35M        35      1322
     EAS1_97:4:274:287:423     seq2      +         35M        35      1332
    EAS54_71:8:105:854:975     seq2      +         35M        35      1354
  EAS139_11:7:50:1229:1313     seq2      +         35M        35      1376
     EAS114_26:7:37:79:581     seq2      +         35M        35      1349
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578       219        35         0
      B7_593:4:106:316:452       259        36         0
    EAS54_65:3:321:311:983       262        35         0
       B7_591:5:42:540:501       259        36         0
    EAS192_3:5:223:142:410       269        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1356        35         0
     EAS1_97:4:274:287:423      1366        35         0
    EAS54_71:8:105:854:975      1388        35         0
  EAS139_11:7:50:1229:1313      1410        35         0
     EAS114_26:7:37:79:581      1383        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> # or
> second(galp)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      -         35M        35       185
      B7_593:4:106:316:452     seq1      -         36M        36       224
    EAS54_65:3:321:311:983     seq1      -         35M        35       228
       B7_591:5:42:540:501     seq1      -         36M        36       224
    EAS192_3:5:223:142:410     seq1      -         35M        35       235
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      +         35M        35      1322
     EAS1_97:4:274:287:423     seq2      +         35M        35      1332
    EAS54_71:8:105:854:975     seq2      +         35M        35      1354
  EAS139_11:7:50:1229:1313     seq2      +         35M        35      1376
     EAS114_26:7:37:79:581     seq2      +         35M        35      1349
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578       219        35         0
      B7_593:4:106:316:452       259        36         0
    EAS54_65:3:321:311:983       262        35         0
       B7_591:5:42:540:501       259        36         0
    EAS192_3:5:223:142:410       269        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1356        35         0
     EAS1_97:4:274:287:423      1366        35         0
    EAS54_71:8:105:854:975      1388        35         0
  EAS139_11:7:50:1229:1313      1410        35         0
     EAS114_26:7:37:79:581      1383        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> 
> strandMode(galp)
[1] 1
> first(galp, real.strand=TRUE)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      +         35M        35        36
      B7_593:4:106:316:452     seq1      +         36M        36        49
    EAS54_65:3:321:311:983     seq1      +         35M        35        51
       B7_591:5:42:540:501     seq1      +         36M        36        60
    EAS192_3:5:223:142:410     seq1      +         35M        35        60
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      -         35M        35      1513
     EAS1_97:4:274:287:423     seq2      -         35M        35      1515
    EAS54_71:8:105:854:975     seq2      -         33M        33      1523
  EAS139_11:7:50:1229:1313     seq2      -         35M        35      1528
     EAS114_26:7:37:79:581     seq2      -         35M        35      1533
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578        70        35         0
      B7_593:4:106:316:452        84        36         0
    EAS54_65:3:321:311:983        85        35         0
       B7_591:5:42:540:501        95        36         0
    EAS192_3:5:223:142:410        94        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1547        35         0
     EAS1_97:4:274:287:423      1549        35         0
    EAS54_71:8:105:854:975      1555        33         0
  EAS139_11:7:50:1229:1313      1562        35         0
     EAS114_26:7:37:79:581      1567        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> last(galp, real.strand=TRUE)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      +         35M        35       185
      B7_593:4:106:316:452     seq1      +         36M        36       224
    EAS54_65:3:321:311:983     seq1      +         35M        35       228
       B7_591:5:42:540:501     seq1      +         36M        36       224
    EAS192_3:5:223:142:410     seq1      +         35M        35       235
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      -         35M        35      1322
     EAS1_97:4:274:287:423     seq2      -         35M        35      1332
    EAS54_71:8:105:854:975     seq2      -         35M        35      1354
  EAS139_11:7:50:1229:1313     seq2      -         35M        35      1376
     EAS114_26:7:37:79:581     seq2      -         35M        35      1349
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578       219        35         0
      B7_593:4:106:316:452       259        36         0
    EAS54_65:3:321:311:983       262        35         0
       B7_591:5:42:540:501       259        36         0
    EAS192_3:5:223:142:410       269        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1356        35         0
     EAS1_97:4:274:287:423      1366        35         0
    EAS54_71:8:105:854:975      1388        35         0
  EAS139_11:7:50:1229:1313      1410        35         0
     EAS114_26:7:37:79:581      1383        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> strand(galp)
factor-Rle of length 1572 with 670 runs
  Lengths: 27  1  2  1  3  1  1  3  2  2  3 ...  1  6  1  2  3  2  1  5  1 51
  Values :  +  -  +  -  +  -  +  -  +  -  + ...  +  -  +  -  +  -  +  -  +  -
Levels(3): + - *
> 
> strandMode(galp) <- 2
> first(galp, real.strand=TRUE)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      -         35M        35        36
      B7_593:4:106:316:452     seq1      -         36M        36        49
    EAS54_65:3:321:311:983     seq1      -         35M        35        51
       B7_591:5:42:540:501     seq1      -         36M        36        60
    EAS192_3:5:223:142:410     seq1      -         35M        35        60
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      +         35M        35      1513
     EAS1_97:4:274:287:423     seq2      +         35M        35      1515
    EAS54_71:8:105:854:975     seq2      +         33M        33      1523
  EAS139_11:7:50:1229:1313     seq2      +         35M        35      1528
     EAS114_26:7:37:79:581     seq2      +         35M        35      1533
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578        70        35         0
      B7_593:4:106:316:452        84        36         0
    EAS54_65:3:321:311:983        85        35         0
       B7_591:5:42:540:501        95        36         0
    EAS192_3:5:223:142:410        94        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1547        35         0
     EAS1_97:4:274:287:423      1549        35         0
    EAS54_71:8:105:854:975      1555        33         0
  EAS139_11:7:50:1229:1313      1562        35         0
     EAS114_26:7:37:79:581      1567        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> last(galp, real.strand=TRUE)
GAlignments object with 1572 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     seq1      -         35M        35       185
      B7_593:4:106:316:452     seq1      -         36M        36       224
    EAS54_65:3:321:311:983     seq1      -         35M        35       228
       B7_591:5:42:540:501     seq1      -         36M        36       224
    EAS192_3:5:223:142:410     seq1      -         35M        35       235
                       ...      ...    ...         ...       ...       ...
  EAS139_11:5:52:1278:1478     seq2      +         35M        35      1322
     EAS1_97:4:274:287:423     seq2      +         35M        35      1332
    EAS54_71:8:105:854:975     seq2      +         35M        35      1354
  EAS139_11:7:50:1229:1313     seq2      +         35M        35      1376
     EAS114_26:7:37:79:581     seq2      +         35M        35      1349
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578       219        35         0
      B7_593:4:106:316:452       259        36         0
    EAS54_65:3:321:311:983       262        35         0
       B7_591:5:42:540:501       259        36         0
    EAS192_3:5:223:142:410       269        35         0
                       ...       ...       ...       ...
  EAS139_11:5:52:1278:1478      1356        35         0
     EAS1_97:4:274:287:423      1366        35         0
    EAS54_71:8:105:854:975      1388        35         0
  EAS139_11:7:50:1229:1313      1410        35         0
     EAS114_26:7:37:79:581      1383        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> strand(galp)
factor-Rle of length 1572 with 670 runs
  Lengths: 27  1  2  1  3  1  1  3  2  2  3 ...  1  6  1  2  3  2  1  5  1 51
  Values :  -  +  -  +  -  +  -  +  -  +  - ...  -  +  -  +  -  +  -  +  -  +
Levels(3): + - *
> 
> seqnames(galp)
factor-Rle of length 1572 with 2 runs
  Lengths:  722  850
  Values : seq1 seq2
Levels(2): seq1 seq2
> 
> head(njunc(galp))
[1] 0 0 0 0 0 0
> table(isProperPair(galp))

TRUE 
1572 
> seqlevels(galp)
[1] "seq1" "seq2"
> 
> ## Rename the reference sequences:
> seqlevels(galp) <- sub("seq", "chr", seqlevels(galp))
> seqlevels(galp)
[1] "chr1" "chr2"
> 
> galp[[1]]
GAlignments object with 2 alignments and 0 metadata columns:
      seqnames strand       cigar    qwidth     start       end     width
         <Rle>  <Rle> <character> <integer> <integer> <integer> <integer>
  [1]     chr1      +         35M        35        36        70        35
  [2]     chr1      -         35M        35       185       219        35
          njunc
      <integer>
  [1]         0
  [2]         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> unlist(galp)
GAlignments object with 3144 alignments and 0 metadata columns:
                           seqnames strand       cigar    qwidth     start
                              <Rle>  <Rle> <character> <integer> <integer>
     EAS54_61:4:143:69:578     chr1      +         35M        35        36
     EAS54_61:4:143:69:578     chr1      -         35M        35       185
      B7_593:4:106:316:452     chr1      +         36M        36        49
      B7_593:4:106:316:452     chr1      -         36M        36       224
    EAS54_65:3:321:311:983     chr1      +         35M        35        51
                       ...      ...    ...         ...       ...       ...
    EAS54_71:8:105:854:975     chr2      +         35M        35      1354
  EAS139_11:7:50:1229:1313     chr2      -         35M        35      1528
  EAS139_11:7:50:1229:1313     chr2      +         35M        35      1376
     EAS114_26:7:37:79:581     chr2      -         35M        35      1533
     EAS114_26:7:37:79:581     chr2      +         35M        35      1349
                                 end     width     njunc
                           <integer> <integer> <integer>
     EAS54_61:4:143:69:578        70        35         0
     EAS54_61:4:143:69:578       219        35         0
      B7_593:4:106:316:452        84        36         0
      B7_593:4:106:316:452       259        36         0
    EAS54_65:3:321:311:983        85        35         0
                       ...       ...       ...       ...
    EAS54_71:8:105:854:975      1388        35         0
  EAS139_11:7:50:1229:1313      1562        35         0
  EAS139_11:7:50:1229:1313      1410        35         0
     EAS114_26:7:37:79:581      1567        35         0
     EAS114_26:7:37:79:581      1383        35         0
  -------
  seqinfo: 2 sequences from an unspecified genome
> 
> grglist(galp)  # a GRangesList object
GRangesList object of length 1572:
$EAS54_61:4:143:69:578 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
         <Rle>  <IRanges>  <Rle>
  [1]     chr1 [185, 219]      -
  [2]     chr1 [ 36,  70]      -

$B7_593:4:106:316:452 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
  [1]     chr1 [224, 259]      -
  [2]     chr1 [ 49,  84]      -

$EAS54_65:3:321:311:983 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
  [1]     chr1 [228, 262]      -
  [2]     chr1 [ 51,  85]      -

...
<1569 more elements>
-------
seqinfo: 2 sequences from an unspecified genome
> 
> strandMode(galp) <- 1
> grglist(galp)
GRangesList object of length 1572:
$EAS54_61:4:143:69:578 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
         <Rle>  <IRanges>  <Rle>
  [1]     chr1 [ 36,  70]      +
  [2]     chr1 [185, 219]      +

$B7_593:4:106:316:452 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
  [1]     chr1 [ 49,  84]      +
  [2]     chr1 [224, 259]      +

$EAS54_65:3:321:311:983 
GRanges object with 2 ranges and 0 metadata columns:
      seqnames     ranges strand
  [1]     chr1 [ 51,  85]      +
  [2]     chr1 [228, 262]      +

...
<1569 more elements>
-------
seqinfo: 2 sequences from an unspecified genome
> 
> stopifnot(identical(unname(elementNROWS(grglist(galp))), njunc(galp) + 2L))
> 
> granges(galp)  # a GRanges object
GRanges object with 1572 ranges and 0 metadata columns:
                           seqnames       ranges strand
                              <Rle>    <IRanges>  <Rle>
     EAS54_61:4:143:69:578     chr1    [36, 219]      +
      B7_593:4:106:316:452     chr1    [49, 259]      +
    EAS54_65:3:321:311:983     chr1    [51, 262]      +
       B7_591:5:42:540:501     chr1    [60, 259]      +
    EAS192_3:5:223:142:410     chr1    [60, 269]      +
                       ...      ...          ...    ...
  EAS139_11:5:52:1278:1478     chr2 [1322, 1547]      -
     EAS1_97:4:274:287:423     chr2 [1332, 1549]      -
    EAS54_71:8:105:854:975     chr2 [1354, 1555]      -
  EAS139_11:7:50:1229:1313     chr2 [1376, 1562]      -
     EAS114_26:7:37:79:581     chr2 [1349, 1567]      -
  -------
  seqinfo: 2 sequences from an unspecified genome
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>