More than 5 years have passed since last update.
UCSC Genome Browser をRのplotで再現する
14
Last updated at Posted at 2013-01-13
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<-65921878to<-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.seqnamesX.startX.endX.widthX.strandX.featureX.groupX.idX.density1chrX65931358659325871230*unknown1CpG:1171読み込んだデータの可視化
各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!
Register as a new user and use Qiita more conveniently
- You get articles that match your needs
- You can efficiently read back useful information
- You can use dark theme