NGSデータ解析まとめ

サカナ研究者の手探りNGS解析(おもに進化生物学)

Haplotype networkを作成するソフトPopART

系統地理などでよく使うハプロタイプネットワークの作成で、PopARTというソフトが使いやすい。

http://popart.otago.ac.nz/documentation.shtml

名前の通り、原色のpopなHPである。

ネットワーク図については、minimum spaning networkやmedian joining network, ほかTCSのstatistical parsimony networkなど、6種類くらいのものに対応している。Input fileはNexus形式で、Sequence名についても文字数の制約などはないようだ。

計算結果のネットワーク図は、画面上で直感的に修正できる。こういうソフトは意外となかったので、かなり便利。HP上でのdocumantationも充実している。

blasrのインストール(cmakeを使用する)

ずいぶん久しぶりの更新(100日振り!)。でも忘れているわけではないのです・・・
Linuxubuntu 16.04 LTS)にblasrをインストールした。いろいろ苦戦したので、方法の簡単なメモを残しておく。
blasrについては、以前も別のコンピュータにインストールしことがあるが、その時といろいろ違っていたので、改めて記録しておく(またすぐに方法が変わる可能性もあるけれど・・・)。

今回は、前回とは異なり、cmakeを使う方法でインストールした。普通にmakeを使う方法は、自分の環境ではうまく行かなかった。

(1) 下準備
最初に、blasrはPython 2.7で動くので、LinuxPythonが3系になっている場合は、2.7にしておく。以下を参照。
Phython3.6をインストールして切り替える - NGSデータ解析まとめ

cmakeをインストールする。以下のサイト参照にした。普通にapt-get installで問題なくインストールできた。
次に、Ninjaをインストールする。以下のコマンドによる。

sudo apt-get update
sudo apt-get install ninja-build

次に、BOOSTをインストール。このサイトを参照に、sudo apt-get install で普通にインストールできる。
さらに、htslibをインストールする。Githubから最新版をダウンロードしてコンパイルする。

(2) HDF5(libhdf5)のインストール
blasrをコンパイルするために、HDF5が必要なのでインストールする。また、Ubuntu 16.04の場合は、このサイトを参照に、libhdf5をインストールでいいと思う(たぶん)。

以下のようにする。

sudo apt-get install libhdf5-dev

(3) blasrのインストール

インストール方法は、以下のサイトを参考にした。
blasr/INSTALL_CMAKE.md at master · PacificBiosciences/blasr · GitHub

まず最初に、適当なpath(ここでは/home/hashi-linux2/)に"blasr_install"ディレクトリを作成し、gitでblasrをダウンロードする。ダウンロードした"blasr"ディレクトリに入る。

mkdir blasr_install
cd blasr_install
git clone https://github.com/PacificBiosciences/blasr.git --recursive
cd blasr

次に、gitでディレクトリ内容のアップデートを行う。

git submodule update --init --remote

ここで、cmake & ninjaによるcompile, buildを行う。

mkdir build && cd build    # build ディレクトリの作成、その中に移動
cmake -GNinja .. && ninja    # compile & build

うまくコンパイルが終了したら、ディレクトリ内に"blasr"のバイナリができる。ここでエラーが出る場合、メッセージにしたがって必要なライブラリをインストールしたりする必要がある。今回、特に引っかかったのは、このサイトにあるタイプのエラーで、このエラーが出た場合、"build"内にできたすべてのファイル、フォルダを削除して再度compileを行うとうまくいった。

ここまでできたら、今度は環境変数LD_LIBRARY_PATHを設定する。以下をターミナルに入力する。また、設定ファイル.profileにも書いておく(LD_LIBRARY_PATHは、それぞれの環境によって異なるので、適宜書き換えて下さい)。

export LD_LIBRARY_PATH="/home/hashi-linux2/blasr_install/blasr/libcpp/alignment:/home/hashi-linux2/blasr_install/blasr/libcpp/hdf:/home/hashi-linux2/blasr_install/blasr/libcpp/pbdata:/usr/lib/x86_64-linux-gnu/hdf5/serial/:/usr/local/lib/"
# メモ:最後の/usr/local/lib/は"libhts.so.2"の位置(環境によって場所は異なる)

動作確認のため、以下のコマンドを入力する。

./blasr --version

これでusageのメッセージが出れば、compileがうまくいっているので、以下のコマンドでインストールする。

sudo ninja install

これで終了! やれやれ。

Python3.6をインストールして切り替える

NGS解析のとき、使用するLinux (Ubuntu 16.04)でPyhton 3.6とPyhton 2.7を切り替えることがあるので、その方法についてのメモ。
pyenvを使用する。

以下のサイトを参考にした。pyenvのインストールから、python3.6.0のインストールまではこちら。

Python3.6.0をUbuntu16.04に導入する。

使用するPythonのバージョンを3.6.0に切り替えるには、以下のようにする。

pyenv global 3.6.0

Pythonのバージョンを確認するには、以下。

python -V

ここでバージョンが切り替わっていれば、"Python3.6.0"と表示される。

再び、バージョンを2.7(システムのデフォルト設定)に戻すには、以下のようにする。

pyenv global system

使用できるPythonのリストは、以下で表示される。

pyenv versions

MacでフォーマットしたHDDをLinux (Ubuntu12.04)でマウントする(メモ)

HDDをフォーマットするときにしばしば迷うのでメモ。

MacでのHDDのフォーマット(MacOS拡張(ジャーナリング))は、LinuxUbuntu)で普通に読み書き可能でマウントできるらしい。以下のURLによると、「MacOS拡張(ジャーナリング)」は、HFS+というフォーマット形式と同じとのこと。なるほど〜。

qiita.com

GATKを用いた多型サイトの抽出:-L optionについて

参考URL: https://www.broadinstitute.org/gatk/guide/article?id=4133


GATK HaplotypeCallerを使用して、リファレンスにマッピングされた特定の個人のヒトゲノムからSNPやINDELなどの多型を抽出するとき、-L optionを使用することで、計算量を節約できます。

たとえば、ヒトゲノムの20番染色体のみから多型を抽出する場合、-L optionの引数を以下のように指定します。

# 変数定義
sample=C8W82ANXX_PG1144_01A15_H1
echo ${sample}

# GATK HaplotypeCaller実行
java -Xmx4g -jar GenomeAnalysisTK.jar \
	-rf BadCigar -rf FailsVendorQualityCheck -rf MappingQualityUnavailable \
	-T HaplotypeCaller -R human_g1k_v37_decoy.fasta \
	-I ${sample}.fixmate.dedup.realign.recal.sort.addRG.bam \
	--dbsnp dbsnp_138.b37.vcf \
	--emitRefConfidence GVCF \
	--variant_index_type LINEAR \
	--variant_index_parameter 128000 \
        -L 20 \        # -L 調べたい染色体番号(scaffold, contig)を引数に指定
	-o ${sample}_raw_variants.20.g.vcf

また、たとえばexomeのみで多型を調べたい場合、Broad Instituteが提供している"Broad.human.exome.b37.interval_list"などを使用することで、exon領域のみの多型を抽出できます("Broad.human.exome.b37.interval_list"についてはこちらを参照)。

# 変数定義
sample=C8W82ANXX_PG1144_01A15_H1
echo ${sample}

# GATK HaplotypeCaller実行
java -Xmx4g -jar GenomeAnalysisTK.jar \
	-rf BadCigar -rf FailsVendorQualityCheck -rf MappingQualityUnavailable \
	-T HaplotypeCaller -R human_g1k_v37_decoy.fasta \
	-I ${sample}.fixmate.dedup.realign.recal.sort.addRG.bam \
	--dbsnp dbsnp_138.b37.vcf \
	--emitRefConfidence GVCF \
	--variant_index_type LINEAR \
	--variant_index_parameter 128000 \
        -L Broad.human.exome.b37.interval_list \        # -L Broad.human.exome.b37.interval_listを引数に指定
	-o ${sample}_raw_variants.exome.g.vcf

もちろん、全ゲノム領域を対象にして多型を抽出したい場合は、-L optionなしで解析します。