GRanges object with a DataFrame as elementMetadata.
...
Parameters passed to control lines in top plot.
stat.y
integer (variable position starting in DataFrame of data, start from
1) or strings (variable names) which indicate the column names.
stat.ylab
y label for stat track(the top track).
stat.assay
default 1L, element of assays.
sig
a character of element meta data column of logical value, indicates
which row is signficant. and will be shown in link lines and rectangle.
sig.col
colors for significant, valid when you specify "sig" argument, the
first color indicates FALSE, non-significant, the second
color indicates TRUE.
stat.coord.trans
transformation used for top plot.
annotation
A list of ggplot object.
width.ratio
Control the segment length of statistic layer.
theme.stat
top plot theme.
theme.align
alignment themes.
linetype
linetype
heights
Heights of each track.
Details
Inspired by some graphics produced in some other packages, for example
in package DEXseq, the author provides graphics with gene
models and linked to an even spaced statistics summary. This is useful
because we always plot everything along the genomic coordinates, but
genomic features like exons are not evenly distributed, so we could
actually treat the statistics associated with exons like categorical
data, and show them as "Paralell Coordinates Plots". This is one
special layout which represent the data in a nice manner and also keep
the genomic structure information. With abliity of tracks,
it's possible to generate such type of a graphic along with other
annotations.
The data we want is a normal GRanges object, and make sure
the intervals are not overlaped with each other(currently), and you
may have multiple columns which store the statistics for multiple
samples, then we produce the graphic we introduced above and users
could pass other annotation track in the function which will be shown
below the main linked track.
The reason you need to pass annotation into the function instead of
binding them by tracks later is because binding manually
with annotation tracks is tricky and this function doesn't return a
ggplot object.
Value
return a frame grob; side-effect (plotting) if plot=T.
Author(s)
Tengfei Yin
Examples
library(TxDb.Hsapiens.UCSC.hg19.knownGene)
library(ggbio)
data(genesymbol, package = "biovizBase")
txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene
model <- exonsBy(txdb, by = "tx")
model17 <- subsetByOverlaps(model, genesymbol["RBM17"])
exons <- exons(txdb)
exon17 <- subsetByOverlaps(exons, genesymbol["RBM17"])
## reduce to make sure there is no overlap
## just for example
exon.new <- reduce(exon17)
## suppose
values(exon.new)$sample1 <- rnorm(length(exon.new), 10, 3)
values(exon.new)$sample2 <- rnorm(length(exon.new), 10, 10)
values(exon.new)$score <- rnorm(length(exon.new))
values(exon.new)$significant <- sample(c(TRUE,FALSE), size = length(exon.new),replace = TRUE)
plotRangesLinkedToData(exon.new, stat.y = c("sample1", "sample2"))
plotRangesLinkedToData(exon.new, stat.y = 1:2)
plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4)
plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4,
sig = "significant")
plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4,
sig = "significant", sig.col = c("gray90","red"))
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(ggbio)
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: ggplot2
Need specific help about ggbio? try mailing
the maintainer or visit http://tengfei.github.com/ggbio/
Attaching package: 'ggbio'
The following objects are masked from 'package:ggplot2':
geom_bar, geom_rect, geom_segment, ggsave, stat_bin, stat_identity,
xlim
Warning message:
replacing previous import 'ggplot2::Position' by 'BiocGenerics::Position' when loading 'ggbio'
> png(filename="/home/ddbj/snapshot/RGM3/R_BC/result/ggbio/plotRangesLinkedToData.Rd_%03d_medium.png", width=480, height=480)
> ### Name: plotRangesLinkedToData
> ### Title: Plot Ranges Linked with Data
> ### Aliases: plotRangesLinkedToData
> ### plotRangesLinkedToData,RangedSummarizedExperiment-method
> ### plotRangesLinkedToData,GenomicRangesORGRangesList-method
>
> ### ** Examples
>
> library(TxDb.Hsapiens.UCSC.hg19.knownGene)
Loading required package: GenomicFeatures
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")'.
> library(ggbio)
> data(genesymbol, package = "biovizBase")
> txdb <- TxDb.Hsapiens.UCSC.hg19.knownGene
> model <- exonsBy(txdb, by = "tx")
> model17 <- subsetByOverlaps(model, genesymbol["RBM17"])
> exons <- exons(txdb)
> exon17 <- subsetByOverlaps(exons, genesymbol["RBM17"])
> ## reduce to make sure there is no overlap
> ## just for example
> exon.new <- reduce(exon17)
> ## suppose
> values(exon.new)$sample1 <- rnorm(length(exon.new), 10, 3)
> values(exon.new)$sample2 <- rnorm(length(exon.new), 10, 10)
> values(exon.new)$score <- rnorm(length(exon.new))
> values(exon.new)$significant <- sample(c(TRUE,FALSE), size = length(exon.new),replace = TRUE)
>
> plotRangesLinkedToData(exon.new, stat.y = c("sample1", "sample2"))
Scale for 'y' is already present. Adding another scale for 'y', which will
replace the existing scale.
> plotRangesLinkedToData(exon.new, stat.y = 1:2)
Scale for 'y' is already present. Adding another scale for 'y', which will
replace the existing scale.
> plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4)
Scale for 'y' is already present. Adding another scale for 'y', which will
replace the existing scale.
> plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4,
+ sig = "significant")
Scale for 'y' is already present. Adding another scale for 'y', which will
replace the existing scale.
> plotRangesLinkedToData(exon.new, stat.y = 1:2, size = 3, linetype = 4,
+ sig = "significant", sig.col = c("gray90","red"))
Scale for 'y' is already present. Adding another scale for 'y', which will
replace the existing scale.
>
>
>
>
>
> dev.off()
null device
1
>