macでインフォマティクス

macでインフォマティクス

HTS (NGS) 関連のインフォマティクス情報についてまとめています。

染色体イデオグラムや染色体間のシンテニープロットを描画する RIdeogram

CRAN 2020 PeerJ karyogram synteny_block 結果の視覚化 (visualization)

 

　近年、シーケンシング技術の発展、特にPacific Biosciences（Eidら、2009）およびOxford Nanopore Technologies（Laverら、2015）を含む第三世代シーケンシングの急速な進歩、BioNanoゲノムマッピング（Caoら、2014）およびハイスループットクロマチンコンフォメーションキャプチャーシーケンシング（Dekkerら、2002）に伴い、より多くの種がそのゲノムを染色体レベルでシーケンシングまたは更新している（Jiao & Schneeberger、2017; Phillippy、2017）。染色体レベルでのゲノムの完成後、いくつかのゲノム特性の概要は、ヒマワリゲノム全体の遺伝子およびトランスポゾンの分布など、種のゲノムをよりよく理解するのに役立つ（Badouin et al)。

　カリオタイプとしても知られるイデオグラムは、真核細胞の核内の染色体の表現型的外観として定義され、2007年に最初のウェブサーバであるIdiographicaがオンラインになって以来、ゲノム全体のデータを可視化するために広く利用されてきた（Kin & Ono, 2007）。環状ゲノム可視化のためのツールは数多く開発されており、Perl言語ベースのツールCircosが最もよく使われている(Krzywinski et al., 2009; Parveen, Khurana & Kumar, 2019)。対照的に、イデオグラム上の全ゲノム情報の非環状プロットのための代替手段は多くない。GenomeGraphs（Durinckら、2009）、gbio（Yin, Cook & Lawrence、2012）、IdeoViz（Pai & Ren、2014）、chromPlot（Orostica & Verdugo、2016）およびchromDraw（Janecka & Lysak、2016）のようなRパッケージ、そしてJavaScriptライブラリのIdeogram.js（Weitzら、2017）およびkaryotypeSV（KaryotypeSV）のような非環状ゲノム可視化ライブラリが開発されてきたが、それらはいくつかの種やデータ可視化タイプに限定されているか、十分なカスタマイズ性に欠けているかのいずれかである。最近では、容量を強化した2つのRパッケージ、karyoploteR (Gel & Serra, 2017)とchromoMap (Anand, 2019)が開発されている。

　しかし、これらの非環状プロットのすべてがCircosのように達成できない機能の1つは、イデオグラム上のベジェ曲線を用いて2つ以上の種の間の関係を可視化することである。この機能は非常に便利で、特に全ゲノム重複の可視化において、ゲノム全体の関係をより直感的に解釈することができる。実際、Circosは通常、ベジェ曲線を用いた種間および種内ゲノム比較の両方において、シンテニーブロックを示すために使用されている（Huら、2019; Wangら、2019）。このように、非環状ゲノムの可視化のためのRパッケージが不足しており、イデオグラム上のベジェ曲線を用いて2つ以上の種間のゲノムワイドな関係を可視化することを可能にしている。

　Scalable Vector Graphics (SVG)は、2次元グラフィックスアプリケーションや画像を記述するための言語である。SVG グラフィックスは eXtensible Markup Language (XML) テキストファイルで定義されているので、どんなテキストエディタや描画ソフトウェアを使っても SVG グラフィックスを簡単に作成・編集できる。多くのRグラフィックスパッケージは、従来のグラフィックスシステムとグリッドグラフィックスシステムの2つのグラフィックスシステム上に構築されている。この論文では、SVG言語に基づいて、全ゲノム情報をイディオグラム上に可視化・マッピングすることができる、種の制限のない高品質なイディオグラムを描画するためのRパッケージ（RIdeogram）を開発した。また、RIdeogramは、イデオグラム上の対称ブロックを結ぶベジェ曲線を用いて、ゲノムの対称性を表示することも可能である。

 

usage

https://cran.r-project.org/web/packages/RIdeogram/RIdeogram.pdf

https://cran.r-project.org/web/packages/RIdeogram/vignettes/RIdeogram.html

 

インストール

Github

#devtools(無い人だけ)
install.packages("devtools")
#devtoolsのロード
library(devtools)
devtools::install_github('TickingClock1992/RIdeogram')

#またはCRANからインストール
install.packages("RIdeogram")

 

 

実行方法

 

ライブラリのロード（*1）。

require(RIdeogram)

 

B、染色体イデオグラムプロット

１、karyotypeデータの読み込み。

#デモデータのロード
data(human_karyotype, package="RIdeogram")

#実際のデータはファイルから読み込む
human_karyotype <- read.table("karyotype.txt", sep = "\t", header = T, stringsAsFactors = F)

読み込めるのはファイルは５カラム構成のファイル（上ではタブ区切りで読み込んでいる）で、最初の列：染色体番号、２列目と３列目：染色体のスタート（０）とエンドのポジション、４列目と５列目：セントロメアのスタートとエンドのポジション、のフォーマットのカリオタイプファイル。自分の種のセントロメア情報がわからない場合、karyotypeファイルには3つの列だけ記載すればよい。

テストデータを見てみる。

> head(human_karyotype)

f:id:kazumaxneo:20210130165428p:plain

 

２、遺伝子密度ファイル（ヒートマップデータ）の読み込み。

#デモデータのロード
data(gene_density, package="RIdeogram")

#実際のデータはファイルから読み込む
gene_density <- read.table("gene_density.txt", sep = "\t", header = T, stringsAsFactors = F)

gene_densityファイルにはラベル情報が含まれており、４つの列からなる。すなわち、最初の列：染色体番号、２列目と３列目：ある領域のスタート（０）とエンドのポジション、４列目：その領域に存在する遺伝子数、というフォーマットのファイルとなる。

f:id:kazumaxneo:20210130165906p:plain

 

３、表示したいフィーチャー（遺伝子などのトラックラベルデータ）の読み込み。

#デモデータのロード
data(Random_RNAs_500, package="RIdeogram")

#実際のデータはファイルから読み込む
Random_RNAs_500 <- read.table("data_2.txt", sep = "\t", header = T, stringsAsFactors = F)

mydata_intervalファイルにはラベル情報が含まれており、6つの列からなる。すなわち、最初の列；ラベルの種類、2番目の列；ラベルの形状（box、triangle、circle）、3番目の列：染色体番号、4番目と5番目の列：染色体内の対応するラベルの開始位置と終了位置、6番目の列：ラベルの色、というフォーマットのファイルとなる。

f:id:kazumaxneo:20210130174841p:plain

 

＜遺伝子密度ファイルがない場合＞

GFFファイルからヒートマップ情報（遺伝子密度など）を抽出するための関数GFFexも用意されている。そのためには、ゲノムのGFFファイルと、そのGFFと同じフォーマットのkaryotypeファイルを用意し、以下のようなコマンドを実行する。

gene_density <- GFFex(input = "gencode.v32.annotation.gff3.gz", karyotype = "human_karyotype.txt", feature = "gene", window = 1000000)

karyotypeファイルの最初の列の染色体IDは、gffファイルのそれと同じでなければならない。

次に、以下のコードを実行。引数 "feature"（デフォルト値は "gene"）でGFFファイルから抽出したい特徴を選択し、引数 "window"（デフォルト値は "1000000"）でウィンドウサイズを設定する。

ideogram(karyotype, overlaid = NULL, label = NULL, label_type = NULL, synteny = NULL, colorset1, colorset2, width, Lx, Ly, output = "chromosome.svg")
convertSVG(svg, device, width, height, dpi)

 

４、染色体イデオグラム（karyogram）を描画する。

ideogram(karyotype = human_karyotype)
#SVGからPNGにconvertSVG関数を使って変換
convertSVG("chromosome.svg", device = "png") 

#TIFF変換。DPIは600とする。 ("tiff", "pdf", "jpg")
convertSVG("chromosome.svg", device = "tiff", dpi = 600)

chromosome.svg出力

f:id:kazumaxneo:20210130180309p:plain

 

５、遺伝子密度のヒートマップを染色体イデオグラムの上にオーバーレイ表示する。結果が出るまで数十秒かかる。

ideogram(karyotype = human_karyotype, overlaid = gene_density)
#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

chromosome.svg出力

f:id:kazumaxneo:20210130180641p:plain


６、染色体イデオグラムの横にトラックラベルを表示。

ideogram(karyotype = human_karyotype, label = Random_RNAs_500, label_type = "marker")
#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

chromosome.svg出力

f:id:kazumaxneo:20210130180802p:plain

 

７、遺伝子密度のヒートマップとトラックラベル両方。

ideogram(karyotype = human_karyotype, overlaid = gene_density, label = Random_RNAs_500, label_type = "marker")
#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

chromosome.svg出力

f:id:kazumaxneo:20210130181309p:plain


８、ヒートマップの色を変更したい場合は、引数 'colorset1' を変更する (デフォルトの設定は colorset1 = c("#4575b4", "#ffffbf", "#d73027")。colors() で指定された色名、 16 進カラーを使用できる。

ideogram(karyotype = human_karyotype, overlaid = gene_density, label = Random_RNAs_500, label_type = "marker", colorset1 = c("#fc8d59", "#ffffbf", "#91bfdb"))
#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

chromosome.svg出力

f:id:kazumaxneo:20210130181605p:plain

 

９、染色体数が多い/少ない場合、width（170）の幅を変更する。 例えば"width = 100"ぐらいにする。凡例を移動させたい場合は、引数 'Lx' と 'Ly' を変更する(デフォルト値は "160" と "35")。

ideogram(karyotype = human_karyotype, overlaid = gene_density, label = Random_RNAs_500, label_type = "marker", width = 100, Lx = 80, Ly = 25) 

#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

chromosome.svg出力

f:id:kazumaxneo:20210130182309p:plain

 

 

 

B、ゲノムシンテニープロット

2ゲノムまたは3ゲノム間のゲノムシンテニー結果の可視化ができる。

１、まず２つのゲノムのイデオグラム情報を読み込む。

#デモデータのロード
data(karyotype_dual_comparison, package="RIdeogram")

##実際のデータはファイルから読み込む
human_karyotype <- read.table("karyotype_dual_comparison.txt", sep = "\t", header = T, stringsAsFactors = F)

> head(karyotype_dual_comparison)

f:id:kazumaxneo:20210130184554p:plain

自分のデータをインポートしたい場合は、上記のようにread.table関数を使う。イデオグラムと少しフォーマットが異なり、最初の3列は同じで、4列目：イデオグラムのボディカラー、5列目：種の名前、６列目と７列目：種名ラベルの大きさと色、を記載する。この核型ファイルには、2つのゲノム（種A：ブドウ、種B：ポピュラス）の情報が含まれており、種Aが前方にソートされている（先頭６行をheadで見ているので、ゲノムの１つ目の方のGrapeだけが見えている）。

$speciesを見て染色体数を確認。

table(karyotype_dual_comparison$species)

f:id:kazumaxneo:20210130184829p:plain

 

２、次にシンテニー座標を記載したファイルを読み込む。

#デモデータのロード
data(synteny_dual_comparison, package="RIdeogram")
head(synteny_dual_comparison)

f:id:kazumaxneo:20210130184901p:plain

デュアルゲノムシンテニーファイルでは、最初の3列：種A(ブドウ)のシンテニー位置情報、次の3列：対応するシンテニーブロックの種B(ポピュラス)のシンテニー位置情報、最後の1列：対応するシンテニーブロックを結ぶベジエ曲線の色、のフォーマットのファイルを読み込む。色のついた線をできるだけ最後にソートする。

 

３、視覚化する。

ideogram(karyotype = karyotype_dual_comparison, synteny = synteny_dual_comparison)

#SVGからPNGに変換
convertSVG("chromosome.svg", device = "png")

f:id:kazumaxneo:20210130185004p:plain

 
 3つのゲノムの場合（マニュアル参照）。

ideogram(karyotype = karyotype_ternary_comparison, synteny = synteny_ternary_comparison)
convertSVG("chromosome.svg", device = "png") 

f:id:kazumaxneo:20210130185317p:plain

 


 

convertSVG以外にも、svg2tiff、svg2pdf、svg2jpg、svg2pngが利用できる。

#SVGをPDF変換。DPIは300とする。 
svg2pdf("chromosome.svg", dpi = 300)

  

他にも染色体イデオグラムの隣にラインやポリゴンを表示させる機能などがあり、これは染色体イデオグラムに沿って何かの情報を表示させる時に利用できます。マニュアルを確認して下さい。

引用

RIdeogram: drawing SVG graphics to visualize and map genome-wide data on the idiograms

Zhaodong Hao, Dekang Lv, Ying Ge, Jisen Shi, Dolf Weijers, Guangchuang Yu​, Jinhui Chen​

PeerJ, January 20, 2020

 

参考