Supported by Dr. Osamu Ogasawara and  providing  . |
|
Last data update: 2014.03.03 |
Extract the edges and their ranges from a SplicingGraphs objectDescription
Alternatively, Usage
sgedgesByGene(x, with.exon.mcols=FALSE, with.hits.mcols=FALSE,
keep.dup.edges=FALSE)
sgedgesByTranscript(x, with.exon.mcols=FALSE, with.hits.mcols=FALSE)
## S4 method for signature 'SplicingGraphs'
intronsByTranscript(x)
Arguments
ValueA GRangesList object named with the gene ids and
where the items are grouped by gene (for The items that are being grouped are the splicing graph edges of type
exon and intron (no artificial edges) for When the grouping is by transcript (i.e. for About duplicated edges: A given edge can typically be shared by more than
1 transcript within the same gene, therefore Author(s)H. Pages See AlsoThis man page is part of the SplicingGraphs package.
Please see Examples
## ---------------------------------------------------------------------
## 1. Make SplicingGraphs object 'sg' from toy gene model (see
## '?SplicingGraphs')
## ---------------------------------------------------------------------
example(SplicingGraphs)
sg
## 'sg' has 1 element per gene and 'names(sg)' gives the gene ids.
names(sg)
## ---------------------------------------------------------------------
## 2. sgedgesByGene()
## ---------------------------------------------------------------------
edges_by_gene <- sgedgesByGene(sg)
edges_by_gene
## 'edges_by_gene' has the length and names of 'sg', that is, the names
## on it are the gene ids and are guaranteed to be unique.
## Extract the edges and their ranges for a given gene:
edges_by_gene[["geneB"]]
## Note that edge with global edge id "geneB:3,4" is an intron that
## belongs to transcripts B1 and B2.
edges_by_gene0 <- sgedgesByGene(sg, keep.dup.edges=TRUE)
edges_by_gene0[["geneB"]]
## Note that edge "geneB:3,4" now shows up twice, once for transcript
## B1, and once for transcript B2.
## Keep the "exon metadata columns":
sgedgesByGene(sg, with.exon.mcols=TRUE)
## Note that those cols are set to NA for intronic edges.
## ---------------------------------------------------------------------
## 3. sgedgesByTranscript()
## ---------------------------------------------------------------------
edges_by_tx <- sgedgesByTranscript(sg)
edges_by_tx
## 'edges_by_tx' is typically longer than 'sg'.
## IMPORTANT NOTE: One caveat here is that the names on 'edges_by_tx'
## are the gene ids, not the transcript ids, and thus are typically NOT
## unique!
## Select elements of a given gene:
edges_by_tx["geneB"] # not a good idea
edges_by_tx[names(edges_by_tx) %in% "geneB"] # much better :-)
## Note that edge with global edge id "geneB:3,4" is an intron that
## belongs to transcripts B1 and B2.
## Keep the "exon metadata columns":
sgedgesByTranscript(sg, with.exon.mcols=TRUE)
## Note that those cols are set to NA for intronic edges.
## ---------------------------------------------------------------------
## 4. intronsByTranscript()
## ---------------------------------------------------------------------
in_by_tx <- intronsByTranscript(sg)
in_by_tx
## 'in_by_tx' has the length and names of 'edges_by_tx'. The same
## recommendation applies for selecting elements of a given set of
## genes:
in_by_tx[c("geneB", "geneD")] # not a good idea
in_by_tx[names(in_by_tx) %in% c("geneB", "geneD")] # much better :-)
## ---------------------------------------------------------------------
## 5. Comparing the outputs of unlist(), intronsByTranscript(), and
## sgedgesByTranscript()
## ---------------------------------------------------------------------
ex_by_tx <- unlist(sg)
in_by_tx <- intronsByTranscript(sg)
edges_by_tx <- sgedgesByTranscript(sg)
## A sanity check:
stopifnot(identical(elementNROWS(in_by_tx) + 1L,
elementNROWS(ex_by_tx)))
## 'edges_by_tx' combines 'ex_by_tx' and 'in_by_tx' in a single
## GRangesList object. Sanity check:
stopifnot(identical(elementNROWS(edges_by_tx),
elementNROWS(ex_by_tx) + elementNROWS(in_by_tx)))
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(SplicingGraphs)
Loading required package: GenomicFeatures
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: AnnotationDbi
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: GenomicAlignments
Loading required package: SummarizedExperiment
Loading required package: Biostrings
Loading required package: XVector
Loading required package: Rsamtools
Loading required package: Rgraphviz
Loading required package: graph
Attaching package: 'graph'
The following object is masked from 'package:Biostrings':
complement
Loading required package: grid
Attaching package: 'Rgraphviz'
The following objects are masked from 'package:IRanges':
from, to
The following objects are masked from 'package:S4Vectors':
from, to
Warning messages:
1: replacing previous import 'IRanges::from' by 'Rgraphviz::from' when loading 'SplicingGraphs'
2: replacing previous import 'IRanges::to' by 'Rgraphviz::to' when loading 'SplicingGraphs'
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/SplicingGraphs/sgedgesByGene-methods.Rd_%03d_medium.png", width=480, height=480)
> ### Name: sgedgesByGene-methods
> ### Title: Extract the edges and their ranges from a SplicingGraphs object
> ### Aliases: sgedgesByGene-methods
> ### intronsByTranscript,SplicingGraphs-method sgedgesByTranscript
> ### sgedgesByTranscript,SplicingGraphs-method sgedgesByGene
> ### sgedgesByGene,SplicingGraphs-method
>
> ### ** Examples
>
> ## ---------------------------------------------------------------------
> ## 1. Make SplicingGraphs object 'sg' from toy gene model (see
> ## '?SplicingGraphs')
> ## ---------------------------------------------------------------------
> example(SplicingGraphs)
SplcnG> ## ---------------------------------------------------------------------
SplcnG> ## 1. Load a toy gene model as a TxDb object
SplcnG> ## ---------------------------------------------------------------------
SplcnG>
SplcnG> library(GenomicFeatures)
SplcnG> suppressWarnings(
SplcnG+ toy_genes_txdb <- makeTxDbFromGFF(toy_genes_gff())
SplcnG+ )
Import genomic features from the file as a GRanges object ... OK
Prepare the 'metadata' data frame ... OK
Make the TxDb object ... OK
SplcnG> ## ---------------------------------------------------------------------
SplcnG> ## 2. Compute all the splicing graphs (1 graph per gene) and return them
SplcnG> ## in a SplicingGraphs object
SplcnG> ## ---------------------------------------------------------------------
SplcnG>
SplcnG> ## Extract the exons grouped by transcript:
SplcnG> ex_by_tx <- exonsBy(toy_genes_txdb, by="tx", use.names=TRUE)
SplcnG> ## Extract the transcripts grouped by gene:
SplcnG> tx_by_gn <- transcriptsBy(toy_genes_txdb, by="gene")
SplcnG> sg <- SplicingGraphs(ex_by_tx, tx_by_gn)
SplcnG> sg
SplicingGraphs object with 5 gene(s) and 13 transcript(s)
SplcnG> ## Alternatively 'sg' can be constructed directly from the TxDb
SplcnG> ## object:
SplcnG> sg2 <- SplicingGraphs(toy_genes_txdb) # same as 'sg'
SplcnG> sg2
SplicingGraphs object with 5 gene(s) and 13 transcript(s)
SplcnG> ## Note that because SplicingGraphs objects have a slot that is an
SplcnG> ## environment (for caching the bubbles), they cannot be compared with
SplcnG> ## 'identical()' (will always return FALSE). 'all.equal()' should be
SplcnG> ## used instead:
SplcnG> stopifnot(isTRUE(all.equal(sg2, sg)))
SplcnG> ## 'sg' has 1 element per gene and 'names(sg)' gives the gene ids:
SplcnG> length(sg)
[1] 5
SplcnG> names(sg)
[1] "geneA" "geneB" "geneC" "geneD" "geneE"
SplcnG> ## ---------------------------------------------------------------------
SplcnG> ## 3. Basic manipulation of a SplicingGraphs object
SplcnG> ## ---------------------------------------------------------------------
SplcnG>
SplcnG> ## Basic accessors:
SplcnG> seqnames(sg)
geneA geneB geneC geneD geneE
chrX chrX chrX chrX chrX
Levels: chrX
SplcnG> strand(sg)
geneA geneB geneC geneD geneE
+ - + + +
Levels: + - *
SplcnG> seqinfo(sg)
Seqinfo object with 1 sequence from an unspecified genome; no seqlengths:
seqnames seqlengths isCircular genome
chrX NA NA <NA>
SplcnG> ## Number of transcripts per gene:
SplcnG> elementNROWS(sg)
geneA geneB geneC geneD geneE
2 2 3 4 2
SplcnG> ## The transcripts of a given gene can be extracted with [[. The result
SplcnG> ## is an *unnamed* GRangesList object containing the exons grouped by
SplcnG> ## transcript:
SplcnG> sg[["geneD"]]
GRangesList object of length 4:
[[1]]
GRanges object with 2 ranges and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
<Rle> <IRanges> <Rle> | <integer> <character> <integer> <integer>
[1] chrX [601, 630] + | 10 Dx2 1 1
[2] chrX [666, 675] + | 12 Dx4 2 5
end_SSid
<integer>
[1] 3
[2] 6
[[2]]
GRanges object with 2 ranges and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
[1] chrX [601, 620] + | 9 Dx1 1 1
[2] chrX [651, 700] + | 11 Dx3 2 4
end_SSid
[1] 2
[2] 8
[[3]]
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
[1] chrX [601, 620] + | 9 Dx1 1 1
[2] chrX [666, 675] + | 12 Dx4 2 5
[3] chrX [691, 700] + | 13 Dx5 3 7
end_SSid
[1] 2
[2] 6
[3] 8
...
<1 more element>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
SplcnG> ## See '?plotTranscripts' for how to plot those transcripts.
SplcnG>
SplcnG> ## The transcripts of all the genes can be extracted with unlist(). The
SplcnG> ## result is a *named* GRangesList object containing the exons grouped
SplcnG> ## by transcript. The names on the object are the gene ids:
SplcnG> ex_by_tx <- unlist(sg)
SplcnG> ex_by_tx
GRangesList object of length 13:
$geneA
GRanges object with 1 range and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
<Rle> <IRanges> <Rle> | <integer> <character> <integer> <integer>
[1] chrX [11, 50] + | 2 Ax2 1 1
end_SSid
<integer>
[1] 3
$geneA
GRanges object with 2 ranges and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
[1] chrX [11, 40] + | 1 Ax1 1 1
[2] chrX [71, 100] + | 3 Ax3 2 4
end_SSid
[1] 2
[2] 5
$geneB
GRanges object with 2 ranges and 5 metadata columns:
seqnames ranges strand | exon_id exon_name exon_rank start_SSid
[1] chrX [251, 300] - | 23 Bx1 1 3
[2] chrX [201, 230] - | 20 Bx2 2 6
end_SSid
[1] 1
[2] 4
...
<10 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> sg
SplicingGraphs object with 5 gene(s) and 13 transcript(s)
>
> ## 'sg' has 1 element per gene and 'names(sg)' gives the gene ids.
> names(sg)
[1] "geneA" "geneB" "geneC" "geneD" "geneE"
>
> ## ---------------------------------------------------------------------
> ## 2. sgedgesByGene()
> ## ---------------------------------------------------------------------
> edges_by_gene <- sgedgesByGene(sg)
> edges_by_gene
GRangesList object of length 5:
$geneA
GRanges object with 4 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [11, 50] + | 1 3 geneA:1,3 ex
[2] chrX [11, 40] + | 1 2 geneA:1,2 ex
[3] chrX [41, 70] + | 2 4 geneA:2,4 in
[4] chrX [71, 100] + | 4 5 geneA:4,5 ex
tx_id
<CharacterList>
[1] A1
[2] A2
[3] A2
[4] A2
$geneB
GRanges object with 5 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex B1
[2] chrX [231, 250] - | 3 4 geneB:3,4 in B1,B2
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex B1
[4] chrX [251, 270] - | 2 3 geneB:2,3 ex B2
[5] chrX [216, 230] - | 4 5 geneB:4,5 ex B2
$geneC
GRanges object with 9 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id
[1] chrX [401, 415] + | 1 2 geneC:1,2 ex C1,C2
[2] chrX [416, 465] + | 2 7 geneC:2,7 in C1
[3] chrX [466, 480] + | 7 8 geneC:7,8 ex C1
[4] chrX [416, 440] + | 2 5 geneC:2,5 in C2
[5] chrX [441, 455] + | 5 6 geneC:5,6 ex C2
[6] chrX [456, 480] + | 6 9 geneC:6,9 in C2
[7] chrX [481, 500] + | 9 10 geneC:9,10 ex C2,C3
[8] chrX [421, 430] + | 3 4 geneC:3,4 ex C3
[9] chrX [431, 480] + | 4 9 geneC:4,9 in C3
...
<2 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## 'edges_by_gene' has the length and names of 'sg', that is, the names
> ## on it are the gene ids and are guaranteed to be unique.
>
> ## Extract the edges and their ranges for a given gene:
> edges_by_gene[["geneB"]]
GRanges object with 5 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex
[2] chrX [231, 250] - | 3 4 geneB:3,4 in
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex
[4] chrX [251, 270] - | 2 3 geneB:2,3 ex
[5] chrX [216, 230] - | 4 5 geneB:4,5 ex
tx_id
<CharacterList>
[1] B1
[2] B1,B2
[3] B1
[4] B2
[5] B2
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## Note that edge with global edge id "geneB:3,4" is an intron that
> ## belongs to transcripts B1 and B2.
>
> edges_by_gene0 <- sgedgesByGene(sg, keep.dup.edges=TRUE)
> edges_by_gene0[["geneB"]]
GRanges object with 6 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex
[2] chrX [231, 250] - | 3 4 geneB:3,4 in
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex
[4] chrX [251, 270] - | 2 3 geneB:2,3 ex
[5] chrX [231, 250] - | 3 4 geneB:3,4 in
[6] chrX [216, 230] - | 4 5 geneB:4,5 ex
tx_id
<character>
[1] B1
[2] B1
[3] B1
[4] B2
[5] B2
[6] B2
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## Note that edge "geneB:3,4" now shows up twice, once for transcript
> ## B1, and once for transcript B2.
>
> ## Keep the "exon metadata columns":
> sgedgesByGene(sg, with.exon.mcols=TRUE)
GRangesList object of length 5:
$geneA
GRanges object with 4 ranges and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [11, 50] + | 1 3 geneA:1,3 ex
[2] chrX [11, 40] + | 1 2 geneA:1,2 ex
[3] chrX [41, 70] + | 2 4 geneA:2,4 in
[4] chrX [71, 100] + | 4 5 geneA:4,5 ex
tx_id exon_id exon_name exon_rank start_SSid
<CharacterList> <IntegerList> <CharacterList> <IntegerList> <integer>
[1] A1 2 Ax2 1 1
[2] A2 1 Ax1 1 1
[3] A2 NA NA NA <NA>
[4] A2 3 Ax3 2 4
end_SSid
<integer>
[1] 3
[2] 2
[3] <NA>
[4] 5
$geneB
GRanges object with 5 ranges and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id exon_id
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex B1 23
[2] chrX [231, 250] - | 3 4 geneB:3,4 in B1,B2 NA
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex B1 20
[4] chrX [251, 270] - | 2 3 geneB:2,3 ex B2 22
[5] chrX [216, 230] - | 4 5 geneB:4,5 ex B2 21
exon_name exon_rank start_SSid end_SSid
[1] Bx1 1 3 1
[2] NA NA <NA> <NA>
[3] Bx2 2 6 4
[4] Bx3 1 3 2
[5] Bx4 2 5 4
$geneC
GRanges object with 9 ranges and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id exon_id
[1] chrX [401, 415] + | 1 2 geneC:1,2 ex C1,C2 4
[2] chrX [416, 465] + | 2 7 geneC:2,7 in C1 NA
[3] chrX [466, 480] + | 7 8 geneC:7,8 ex C1 7
[4] chrX [416, 440] + | 2 5 geneC:2,5 in C2 NA
[5] chrX [441, 455] + | 5 6 geneC:5,6 ex C2 6
[6] chrX [456, 480] + | 6 9 geneC:6,9 in C2 NA
[7] chrX [481, 500] + | 9 10 geneC:9,10 ex C2,C3 8
[8] chrX [421, 430] + | 3 4 geneC:3,4 ex C3 5
[9] chrX [431, 480] + | 4 9 geneC:4,9 in C3 NA
exon_name exon_rank start_SSid end_SSid
[1] Cx1 1 1 2
[2] NA NA <NA> <NA>
[3] Cx4 2 7 8
[4] NA NA <NA> <NA>
[5] Cx3 2 5 6
[6] NA NA <NA> <NA>
[7] Cx5 3,2 9 10
[8] Cx2 1 3 4
[9] NA NA <NA> <NA>
...
<2 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## Note that those cols are set to NA for intronic edges.
>
> ## ---------------------------------------------------------------------
> ## 3. sgedgesByTranscript()
> ## ---------------------------------------------------------------------
> edges_by_tx <- sgedgesByTranscript(sg)
> edges_by_tx
GRangesList object of length 13:
$geneA
GRanges object with 1 range and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [11, 50] + | 1 3 geneA:1,3 ex
tx_id
<character>
[1] A1
$geneA
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id
[1] chrX [11, 40] + | 1 2 geneA:1,2 ex A2
[2] chrX [41, 70] + | 2 4 geneA:2,4 in A2
[3] chrX [71, 100] + | 4 5 geneA:4,5 ex A2
$geneB
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex B1
[2] chrX [231, 250] - | 3 4 geneB:3,4 in B1
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex B1
...
<10 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
>
> ## 'edges_by_tx' is typically longer than 'sg'.
> ## IMPORTANT NOTE: One caveat here is that the names on 'edges_by_tx'
> ## are the gene ids, not the transcript ids, and thus are typically NOT
> ## unique!
>
> ## Select elements of a given gene:
> edges_by_tx["geneB"] # not a good idea
GRangesList object of length 1:
$geneB
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex
[2] chrX [231, 250] - | 3 4 geneB:3,4 in
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex
tx_id
<character>
[1] B1
[2] B1
[3] B1
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> edges_by_tx[names(edges_by_tx) %in% "geneB"] # much better :-)
GRangesList object of length 2:
$geneB
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex
[2] chrX [231, 250] - | 3 4 geneB:3,4 in
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex
tx_id
<character>
[1] B1
[2] B1
[3] B1
$geneB
GRanges object with 3 ranges and 5 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id
[1] chrX [251, 270] - | 2 3 geneB:2,3 ex B2
[2] chrX [231, 250] - | 3 4 geneB:3,4 in B2
[3] chrX [216, 230] - | 4 5 geneB:4,5 ex B2
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## Note that edge with global edge id "geneB:3,4" is an intron that
> ## belongs to transcripts B1 and B2.
>
> ## Keep the "exon metadata columns":
> sgedgesByTranscript(sg, with.exon.mcols=TRUE)
GRangesList object of length 13:
$geneA
GRanges object with 1 range and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in
<Rle> <IRanges> <Rle> | <character> <character> <character> <factor>
[1] chrX [11, 50] + | 1 3 geneA:1,3 ex
tx_id exon_id exon_name exon_rank start_SSid end_SSid
<character> <integer> <character> <integer> <integer> <integer>
[1] A1 2 Ax2 1 1 3
$geneA
GRanges object with 3 ranges and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id exon_id
[1] chrX [11, 40] + | 1 2 geneA:1,2 ex A2 1
[2] chrX [41, 70] + | 2 4 geneA:2,4 in A2 <NA>
[3] chrX [71, 100] + | 4 5 geneA:4,5 ex A2 3
exon_name exon_rank start_SSid end_SSid
[1] Ax1 1 1 2
[2] <NA> <NA> <NA> <NA>
[3] Ax3 2 4 5
$geneB
GRanges object with 3 ranges and 10 metadata columns:
seqnames ranges strand | from to sgedge_id ex_or_in tx_id exon_id
[1] chrX [251, 300] - | 1 3 geneB:1,3 ex B1 23
[2] chrX [231, 250] - | 3 4 geneB:3,4 in B1 <NA>
[3] chrX [201, 230] - | 4 6 geneB:4,6 ex B1 20
exon_name exon_rank start_SSid end_SSid
[1] Bx1 1 3 1
[2] <NA> <NA> <NA> <NA>
[3] Bx2 2 6 4
...
<10 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> ## Note that those cols are set to NA for intronic edges.
>
> ## ---------------------------------------------------------------------
> ## 4. intronsByTranscript()
> ## ---------------------------------------------------------------------
> in_by_tx <- intronsByTranscript(sg)
> in_by_tx
GRangesList object of length 13:
$geneA
GRanges object with 0 ranges and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
$geneA
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
[1] chrX [41, 70] +
$geneB
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
[1] chrX [231, 250] -
...
<10 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
>
> ## 'in_by_tx' has the length and names of 'edges_by_tx'. The same
> ## recommendation applies for selecting elements of a given set of
> ## genes:
> in_by_tx[c("geneB", "geneD")] # not a good idea
GRangesList object of length 2:
$geneB
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chrX [231, 250] -
$geneD
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
[1] chrX [631, 665] +
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
> in_by_tx[names(in_by_tx) %in% c("geneB", "geneD")] # much better :-)
GRangesList object of length 6:
$geneB
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
<Rle> <IRanges> <Rle>
[1] chrX [231, 250] -
$geneB
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
[1] chrX [231, 250] -
$geneD
GRanges object with 1 range and 0 metadata columns:
seqnames ranges strand
[1] chrX [631, 665] +
...
<3 more elements>
-------
seqinfo: 1 sequence from an unspecified genome; no seqlengths
>
> ## ---------------------------------------------------------------------
> ## 5. Comparing the outputs of unlist(), intronsByTranscript(), and
> ## sgedgesByTranscript()
> ## ---------------------------------------------------------------------
> ex_by_tx <- unlist(sg)
> in_by_tx <- intronsByTranscript(sg)
> edges_by_tx <- sgedgesByTranscript(sg)
>
> ## A sanity check:
> stopifnot(identical(elementNROWS(in_by_tx) + 1L,
+ elementNROWS(ex_by_tx)))
>
> ## 'edges_by_tx' combines 'ex_by_tx' and 'in_by_tx' in a single
> ## GRangesList object. Sanity check:
> stopifnot(identical(elementNROWS(edges_by_tx),
+ elementNROWS(ex_by_tx) + elementNROWS(in_by_tx)))
>
>
>
>
>
> dev.off()
null device
1
>
|
Created & Maintained by Osamu Ogasawara (osamu.ogasawara@gmail.com) and