USCS Genome Browserは、ゲノム中の遺伝子領域やCpG Island、配列保存度や一塩基多型部位などさまざまな情報（Track）を並べて可視化することができるサイトである。上の図はマウスゲノム(NCBI37/mm9)のX染色体65,921,878 - 65,980,988部分を表示している。この図を得るためには、ふつう上記のサイトで
- group: mammal, genome: Mouse, assembly: July 2007 (NCBI37/mm9) を選択
- テキストボックスに chrX:65921878 - 65980988 と入れてsubmit
- 表示したい個々のTrackを選択してrefresh
のような動作をGUI的に行う必要がある。
これらの一連の処理をRのスクリプトを書いて再現出来れば、ある基準によって抽出したゲノム中の領域について、知りたい個々の情報を自動的にplotすることができるのでとても便利な気がする。
...ということをGvizパッケージを用いるとできるようになるのでメモしておく。
準備
インストール
source("http://bioconductor.org/biocLite.R")
biocLite("Gviz")
パッケージの読み込み
library(Gviz)
表示領域の指定
genome <- "mm9"
chr <- "chrX"
from <- 65921878
to <- 65980988
データの読み込み
各Trackの読み込み
# UCSCにアクセスしているため、一つのTrackを読み込むのに30~40秒くらいかかる
knownGenes <- UcscTrack(genome=genome, chromosome=chr, track="knownGene", from=from, to=to,
trackType="GeneRegionTrack", rstarts="exonStarts", rends="exonEnds",
gene="name", symbol="name", transcript="name", strand="strand",
fill="#8282d2", name="UCSC Genes")
refGenes <- UcscTrack(genome=genome, chromosome=chr, track="xenoRefGene", from=from, to=to,
trackType="GeneRegionTrack", rstarts="exonStarts", rends="exonEnds",
gene="name", symbol="name2", transcript="name", strand="strand",
fill="#8282d2", stacking="dense", name="Other RefSeq")
ensGenes <- UcscTrack(genome=genome, chromosome=chr, track="ensGene", from=from, to=to,
trackType="GeneRegionTrack", rstarts="exonStarts", rends="exonEnds",
gene="name", symbol="name2", transcript="name", strand="strand",
fill="#960000", stacking="dense", name="Ensembl Genes")
cpgIslands <- UcscTrack(genome=genome, chromosome=chr, track="cpgIslandExt", from=from, to=to,
trackType="AnnotationTrack", start="chromStart", end="chromEnd", id="name",
shape="box", fill="#006400", name="CpG Islands")
gcContent <- UcscTrack(genome=genome, chromosome=chr, track="GC Percent", table="gc5Base",
from=from, to=to, trackType="DataTrack", start="start", end="end", data="score",
type="hist", window=-1, windowSize=1500, fill.histogram="black", col.histogram="black",
ylim=c(30, 70), name="GC Percent")
conservation <- UcscTrack(genome=genome, chromosome=chr, track="Conservation", table="phyloP30wayPlacental",
from=from, to=to, trackType="DataTrack", start="start", end="end", data="score",
type="hist", window="auto", col.histogram="darkblue", fill.histogram="darkblue",
ylim=c(-3.7, 4), name="Conservation")
snpLocations <- UcscTrack(genome=genome, chromosome=chr, track="snp128", from=from, to=to,
trackType="AnnotationTrack", start="chromStart", end="chromEnd",
id="name", feature="func", strand="strand", shape="box",
stacking="dense", fill="black", name="SNPs")
axTrack <- GenomeAxisTrack()
idxTrack <- IdeogramTrack(genome=genome, chromosome=chr)
各Trackの中身をデータフレーム形式で閲覧したいときはas()関数を使う。
result
> as(cpgIslands,"data.frame")
X.seqnames X.start X.end X.width X.strand X.feature X.group X.id X.density
1 chrX 65931358 65932587 1230 * unknown 1 CpG: 117 1
読み込んだデータの可視化
各Trackが読み込めたら、最後にplotTracks()でplotする。
plot
plotTracks(list(idxTrack, axTrack, knownGenes,refGenes,ensGenes, cpgIslands,
gcContent, conservation, snpLocations), from=from, to=to, showTitle=TRUE)
なお、ここで読み込んだデータは、実はパッケージ内にデモデータとして組み込まれているので、
UCSCに時間をかけてアクセスするのが面倒な場合は以下のコードを実行すれば
同 じ よ う な 効 果 が 得 ら れ る
忙しい人のためのGvizサンプルコード
library(Gviz)
data(ucscItems)
plotTracks(list(idxTrack, axTrack, knownGenes,refGenes,ensGenes, cpgIslands,
gcContent, conservation, snpLocations),
from=from, to=to, showTitle=TRUE)
くらべてみよう
UCSCで書いた図
configure
- image width: 800px
- Base Position: dense
- Chromosome Band: pack
- GC Percent: full
- UCSC Genes: full
- RefSeq Genes: dense (Hide non-coding genes: check)
- Ensembl Genese: full
- CpG Islands: full
- Conservation: full (Data view scaling: auto-scaling to data view)
- SNPs (128): dense
Gvizパッケージで書いた図
enjoy!