Last data update: 2014.03.03

R: Map range coordinates between proteins and genome space
mapToGenome-methodsR Documentation

Map range coordinates between proteins and genome space

Description

Map range coordinates between peptide features along proteins and genome (reference) space.

Usage

## S4 method for signature 'Proteins,GRangesList'
mapToGenome(x, genome, drop.empty.ranges = TRUE, ...)
## S4 method for signature 'Proteins,GRangesList'
pmapToGenome(x, genome, drop.empty.ranges = TRUE, ...)

Arguments

x

Proteins object containing peptides pranges to be mapped.

genome

A GRangesList object used to map between x and the result. The ranges are typically created by the etrid2grl function.

drop.empty.ranges

TRUE (default) or FALSE. Should empty ranges be dropped?

...

Additional parameters passed to inner functions. Currently ignored.

Details

  • mapToGenome maps the pranges(x) to the ranges of genome. Unless x and genome are of length 1, both must be named and items of x are matched to items of genome using their respective names. Names that do not co-occur in x and genome are ignored. If we have

    seqnames(x): "A", "B" and "C"

    and

    names(genome): "C", "A", "a", "z", "A" and "A".

    the names of the output will be

    "A", "A", "A" and "C".

    The output is ordered by (1) seqnames(x) and (2) the order of the elements in genome.

    In case less than length(x) are mapped, as for p["B"] above, a message informs the user.

  • pmapToGenome is the element-wise (aka 'parallel') version of mapToGenome. The i-th pranges(x) is mapped to the i-th range in genome. x and genome must have the same length and do not need to be named (names are ignored).

Value

A named GRangesList object, with names matching names(genome). For pmapToGenome, the return value will have the same length as the inputs.

Author(s)

Laurent Gatto

See Also

  • See ?mapToAlignments in the GenomicAlignments package for mapping coordinates between reads (local) and genome (reference) space using a CIGAR alignment.

  • See ?mapToTranscripts in the GenomicRanges package for mapping coordinates between features in the transcriptome and genome space.

  • The proteinCoding function to remove non-protein coding ranges before mapping peptides to their genomic coordinates.

  • The mapping vignette for examples and visualisations.

See plotAsAnnotationTrack and plotAsAnnotationTrack for more details about the two plotting functions.

Examples


data(p)

grl <- etrid2grl(acols(p)$ENST)
pcgrl <- proteinCoding(grl)

plotAsGeneRegionTrack(grl[[1]],
                      pcgrl[[1]])

mp <- mapToGenome(p[4], pcgrl[4])

plotAsAnnotationTrack(mp[[1]], pcgrl[[4]])

pmapToGenome(p, pcgrl)

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(Pbase)
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: Rcpp
Loading required package: Gviz
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: GenomicRanges
Loading required package: GenomeInfoDb
Loading required package: grid

This is Pbase version 0.12.2

> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/Pbase/coordinate-mapping-methods.Rd_%03d_medium.png", width=480, height=480)
> ### Name: mapToGenome-methods
> ### Title: Map range coordinates between proteins and genome space
> ### Aliases: mapToGenome mapToGenome-methods mapToGenome,ANY,ANY-method
> ###   mapToGenome,Proteins,GRangesList-method pmapToGenome
> ###   pmapToGenome-methods pmapToGenome,ANY,ANY-method
> ###   pmapToGenome,Proteins,GRangesList-method
> ### Keywords: methods utilities
> 
> ### ** Examples
> 
> 
> data(p)
> 
> grl <- etrid2grl(acols(p)$ENST)
> pcgrl <- proteinCoding(grl)
> 
> plotAsGeneRegionTrack(grl[[1]],
+                       pcgrl[[1]])
> 
> mp <- mapToGenome(p[4], pcgrl[4])
> 
> plotAsAnnotationTrack(mp[[1]], pcgrl[[4]])
> 
> pmapToGenome(p, pcgrl)
GRangesList object of length 9:
$ENST00000409195 
GRanges object with 36 ranges and 7 metadata columns:
    seqnames                 ranges strand |             pepseq   accession
       <Rle>              <IRanges>  <Rle> |        <character> <character>
        chr2 [167249619, 167249672]      + | QEITQNKSFFSSVKESQR      A4UGR9
        chr2 [167239915, 167239950]      + |       LPVPKDVYSKQR      A4UGR9
        chr2 [167246964, 167247002]      + |      EQNNDALEKSLRR      A4UGR9
        chr2 [167246487, 167246516]      + |         SLKESSHRWK      A4UGR9
        chr2 [167249256, 167249303]      + |   LKMVPRKQREFSGSDR      A4UGR9
  .      ...                    ...    ... .                ...         ...
        chr2 [166903540, 166903575]      + |       PESGFAEDSAAR      A4UGR9
        chr2 [167246526, 167246579]      + | QPDAIPGDIEKAIECLEK      A4UGR9
        chr2 [166903624, 166903665]      + |     MARYQAAVSRGDCR      A4UGR9
        chr2 [167247636, 167247674]      + |      TNTSTGLKMAMER      A4UGR9
        chr2 [167249619, 167249660]      + |     QEITQNKSFFSSVK      A4UGR9
               gene      transcript      symbol exonJunctions     group
        <character>     <character> <character>     <logical> <integer>
    ENSG00000163092 ENST00000409195       XIRP2         FALSE         1
    ENSG00000163092 ENST00000409195       XIRP2         FALSE         2
    ENSG00000163092 ENST00000409195       XIRP2         FALSE         3
    ENSG00000163092 ENST00000409195       XIRP2         FALSE         4
    ENSG00000163092 ENST00000409195       XIRP2         FALSE         5
  .             ...             ...         ...           ...       ...
    ENSG00000163092 ENST00000409195       XIRP2         FALSE        32
    ENSG00000163092 ENST00000409195       XIRP2         FALSE        33
    ENSG00000163092 ENST00000409195       XIRP2         FALSE        34
    ENSG00000163092 ENST00000409195       XIRP2         FALSE        35
    ENSG00000163092 ENST00000409195       XIRP2         FALSE        36

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