blastn
1.0.0
01. Blastn データベースのダウンロードと更新 [NCBI ツールによる自動化]
02. Blastn データベースのダウンロードと更新 [bash スクリプトによる自動化]
03. Blastn の実行 [ユーザー対話型 bash スクリプト & Nextflow DSL2 スクリプト]
04. Blastn ヒット配列抽出 [ユーザー対話型 bash スクリプト]
05. Blastn の一般的なエラーと解決策
blastn 用の conda 環境を作成する
conda create -n blastn_db
blastn 環境をアクティブ化する
conda activate blastn_db
次のリンクから最新の blast 実行可能ファイルのリンクをコピーします。
https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/
次のように実行可能ファイルをダウンロードします
wget https://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/ncbi-blast-2.15.0+-x64-linux.tar.gz
ダウンロードしたファイルを次のように解凍します
tar -zxvf ncbi-blast-2.15.0+-x64-linux.tar.gz
次のディレクトリに移動します
cd ncbi-blast-2.15.0+/bin
このパスを PATH 環境変数に追加します。私のチュートリアルでこれを行う方法を参照してください。
https://github.com/asadprodhan/About-the-PATH
update_blastdb.pl を ncbi-blast-2.15.0+/bin ディレクトリから blastn データベースをダウンロードするディレクトリにコピーします
cp ./update_blastdb.pl databaseDirectory
次のようにスクリプトを実行します
run ./update_blastdb.pl --decompress nt
ダウンロードが自動的に開始されます
ダウンロードが完了したら、すべての nt ファイルがダウンロードされたことを確認します。これを行うには、https://ftp.ncbi.nlm.nih.gov/blast/db/ とディレクトリの間で nt ファイル番号をクロスチェックします。
ディレクトリにすべての nt ファイルがない場合は、「BLAST データベース エラー: 参照されたエイリアス ファイルにボリュームまたはエイリアス ファイル nt.xxx が見つかりませんでした」というメッセージが表示されます。
次の bash スクリプトを使用して、不足している nt ファイルをダウンロードできます。
すべての nt ファイルがダウンロードされたら、次のように md5 ファイルを削除できます。
rm -r *.md5
すべての nt.??.tar.gz ファイルのリストを含むmetadata.tsv ファイルを準備します。 nt.??.tar.gz ファイルは次の場所にあります。
https://ftp.ncbi.nlm.nih.gov/blast/db/
ファイルmetadata.tsvは次のようになります。
図 1: Blastn データベースの nt ファイル。
metadata.tsvファイルと次のblastn スクリプトを、blastn データベースをダウンロードするディレクトリに配置します。
次のようにファイル形式を確認してください。
file *
すべてのファイルは UNIX 形式、つまり ASCII テキストのみである必要があります。 Windows コンピュータで書き込まれたファイルは Windows 形式、つまり CRLF 行終端文字を含む ASCII テキストになります。次のコマンドを実行して、これらのファイルを UNIX 形式に変換します。
dos2unix *
ファイルが実行可能であることを確認する
ls -l
次のコマンドを実行してファイルを実行可能にします
chmod +x *
#!/bin/bash
#metadata
metadata=./*.tsv
#
Red="$(tput setaf 1)"
Green="$(tput setaf 2)"
Bold=$(tput bold)
reset=`tput sgr0` # turns off all atribute
while IFS=, read -r field1
do
echo "${Red}${Bold}Downloading ${reset}: "${field1}""
echo ""
wget https://ftp.ncbi.nlm.nih.gov/blast/db/"${field1}"
echo "${Green}${Bold}Downloaded ${reset}: ${field1}"
echo ""
echo "${Green}${Bold}Extracting ${reset}: ${field1}"
tar -xvzf "${field1}"
echo "${Green}${Bold}Extracted ${reset}: ${field1}"
echo ""
echo "${Green}${Bold}Deleting zipped file ${reset}: ${field1}"
rm -r "${field1}"
echo "${Green}${Bold}Deleted ${reset}: ${field1}"
echo ""
done < ${metadata}x: tar にファイルを抽出するように指示します。
v: 「v」は「verbose」を表し、解凍の続行中にすべてのファイルをリストします。
z: tar コマンドにファイル (gzip) を解凍/解凍するように指示します。
f: ファイルが作業対象として割り当てられることを tar に伝えます。
pkill -9 wget # 実行中の wget ダウンロードを中止します
「 tar.gz」のようなワイルドカードは「tar」には機能しません。 「 」で指定された tar は、ディレクトリ内の既存の tar ファイルに限定されるだけでなく、abc.tar.gz def.tar.gz ghi.tar などの架空のファイル名 (!) にも拡張されるためです。 .gz または 1.gz、2.gz、3.gz など。これらのファイルは存在しないため、tar はそれらのファイルを見つけることができず、「アーカイブに見つかりません」エラーが発生します。解凍する tar ファイルが複数ある場合、次のループ関数を使用すると、この問題を解決できます。
for file in *.tar.gz; do tar -xvzf "$file"; done
#!/bin/bash -i
# ask for query file
echo Enter your input file name including extension and hit ENTER
read -e F
# ask for an output directory name
echo Enter an output directory name and hit ENTER
read -e outDir
# ask for the blast database path
echo Enter the path to the blast database and hit ENTER
read -e BlastDB
echo ""
# start monitoring run time
SECONDS=0
# make blast results directory
mkdir ${outDir}
# prepare output file name prefix
baseName=$(basename $F .fasta)
# Run blastn with .asn output
echo blastn in progress...
blastn -db ${BlastDB} -num_alignments 1 -num_threads 16 -outfmt 11 -query $PWD/$F > $PWD/${outDir}/${baseName}.asn
# convert output file from asn to xml format
echo converting output file from asn to xml format
blast_formatter -archive $PWD/${outDir}/${baseName}.asn -outfmt 5 > $PWD/${outDir}/${baseName}.xml
# convert output file from asn to tsv format
echo converting output file from asn to tsv format
blast_formatter -archive $PWD/${outDir}/${baseName}.asn -outfmt 0 > $PWD/${outDir}/${baseName}.tsv
# display the compute time
if (( $SECONDS > 3600 )) ; then
let "hours=SECONDS/3600"
let "minutes=(SECONDS%3600)/60"
let "seconds=(SECONDS%3600)%60"
echo "Completed in $hours hour(s), $minutes minute(s) and $seconds second(s)"
elif (( $SECONDS > 60 )) ; then
let "minutes=(SECONDS%3600)/60"
let "seconds=(SECONDS%3600)%60"
echo "Completed in $minutes minute(s) and $seconds second(s)"
else
echo "Completed in $SECONDS seconds"
fiデフォルトの ID パーセンテージは 90% です
デフォルトのクエリ カバレッジ パーセンテージは 0% です
E 値が小さいほど、一致が良くなります。
参照: https://www.metagenomics.wiki/tools/blast/evalue
ビットスコアが高いほど、シーケンスの類似性が高くなります。
このスクリプトにより、次のことが可能になります
デフォルトの blastn パラメータの変更
コンテナを使用してローカル コンピューターとリモート コンピューターの両方で blastn を実行します (これにより、blastn ソフトウェアのインストールと更新の必要がなくなります。ただし、Nextflow と Singularity をインストールする必要があります)。
複数サンプルの blastn 分析の自動化
#!/usr/bin/env nextflow
nextflow.enable.dsl=2
//data_location
params.in = "$PWD/*.fasta"
params.outdir = './results'
params.db = "./blastn_db"
params.evalue='0.05'
params.identity='90'
params.qcov='90'
// blastn
process blastn {
errorStrategy 'ignore'
tag { file }
publishDir "${params.outdir}/blastn", mode:'copy'
input:
path (file)
path db
output:
path "${file.simpleName}_blast.xml"
path "${file.simpleName}_blast.html"
path "${file.simpleName}_blast_sort_withHeader.tsv"
script:
"""
blastn
-query $file -db ${params.db}/nt
-outfmt 11 -out ${file.simpleName}_blast.asn
-evalue ${params.evalue}
-perc_identity ${params.identity}
-qcov_hsp_perc ${params.qcov}
-num_threads ${task.cpus}
blast_formatter
-archive ${file.simpleName}_blast.asn
-outfmt 5 -out ${file.simpleName}_blast.xml
blast_formatter
-archive ${file.simpleName}_blast.asn
-html -out ${file.simpleName}_blast.html
blast_formatter
-archive ${file.simpleName}_blast.asn
-outfmt "6 qaccver saccver pident length evalue bitscore stitle" -out ${file.simpleName}_blast_unsort.tsv
sort -k1,1 -k5,5n -k4,4nr -k6,6nr ${file.simpleName}_blast_unsort.tsv > ${file.simpleName}_blast_sort.tsv
awk 'BEGIN{print "qaccvertsaccvertpidenttlengthtevaluetbitscoretmismatchtgapopentqstarttqendtsstarttsendtstitle"}1' ${file.simpleName}_blast_sort.tsv > ${file.simpleName}_blast_sort_withHeader.tsv
"""
}
workflow {
query_ch = Channel.fromPath(params.in)
db = file( params.db )
blastn (query_ch, db)
}resume = true
process {
withName:'blastn|blastIndex' { container = 'quay.io/biocontainers/blast:2.14.1--pl5321h6f7f691_0' }
}
singularity {
enabled = true
autoMounts = true
//runOptions = '-e TERM=xterm-256color'
envWhitelist = 'TERM'
}nextflow run main.nf --evalue=0.05 --identity='90' --qcov='0' --db="/path/to/blastn_database"
#!/bin/bash -i
#
Red="$(tput setaf 1)"
Green="$(tput setaf 2)"
Bold=$(tput bold)
reset=`tput sgr0` # turns off all atribute
# ask for blastn output file
echo ""
echo ""
echo "${Red}${Bold}Enter blastn output tsv file and hit ENTER ${reset}"
echo ""
read -e F
echo ""
# ask for the key word
echo "${Red}${Bold}Enter filter word (CASE-SENSITIVE) and hit ENTER ${reset}"
echo ""
read -e KeyWord
echo ""
# ask for the blastn query fasta file
echo "${Red}${Bold}Enter blastn query fasta file and hit ENTER ${reset}"
echo ""
read -e Query
echo ""
# prepare output file name prefix
baseName=$(basename $F .tsv)
echo ""
# filtering the selected blastn hits
echo ""
echo "${Green}${Bold}Filtering the blastn hits containing ${reset}: "${KeyWord}""
echo ""
grep ${KeyWord} $F > ${baseName}_${KeyWord}.tsv
# collecting the query IDs from the selected blastn hits
echo "${Green}${Bold}Collecting the query IDs from the selected blastn hits ${reset}: "${KeyWord}""
echo ""
awk '{print $1}' ${baseName}_${KeyWord}.tsv > IDs.txt
# extracting the sequences for the selected blastn hits
echo "${Green}${Bold}Extracting the sequences for the selected blastn hits ${reset}: "${KeyWord}""
echo ""
bioawk -cfastx 'BEGIN{while((getline k <"IDs.txt")>0)i[k]=1}{if(i[$name])print ">"$name"n"$seq}' ${Query} > ${baseName}_${KeyWord}.fasta
echo ""
echo "${Green}${Bold}Done ${reset}: "${KeyWord}""
echo ""
echo ""blastn ヒット シーケンスを抽出する bash スクリプトはユーザー対話型です。入力を要求し、それらを自動的に処理して、予期される fasta シーケンスを含むファイルを生成します。以下のスクリーンショットを参照してください。
図 2: blastn_hits_sequences_extraction_auto_AP スクリプトの仕組み。
図 3: Blastn データベース エラー「エイリアスまたはインデックス ファイルが見つかりません」。
このエラーは、次のようにスクリプトを調整することで解決できる可能性があります。
/path/to/the/blastn/db/nt のように、データベース パスの末尾に「nt」を追加します。
上記の blastn スクリプトのデータベース パスを参照してください。同様に、blastp の場合は「/nr」
Blastn が Nextflow スクリプトの最初または唯一のプロセスである場合。その場合、プロセスはデータベースのパスを使用する可能性があります。そうでない場合は、データベースをファイルとして提供する必要があります。以下の参考文献を参照してください。また、入力チャネルには path(query_sequence) に加えて path(db) が必要です。上記の blastn スクリプトを参照してください。
https://stackoverflow.com/questions/75465741/path-not-being-detected-by-nextflow
図 4: Blastn データベース エラー「有効なバージョン 4 データベースではありません」。
これは Blast バージョンの競合です
conda 環境を作成すると、最新の blast nr データベースを使用できない blast v2.6 が自動的にインストールされます。
最新の blast nr データベースを使用するには、blast v2.15.0 などの日付のないバージョンが必要です。
使用している blastn のバージョンを確認します。
blastn -version
blastn の最新バージョンを更新します。
conda install -c bioconda blast
図 5: Blastn データベース エラー「参照エイリアスに nt.XXX エイリアスが見つかりませんでした」。
blastn データベース ディレクトリにすべての nt ファイルがない場合、「BLAST データベース エラー: 参照されたエイリアス ファイルにボリュームまたはエイリアス ファイル nt.xxx が見つかりませんでした」というエラーが表示されます。
https://ftp.ncbi.nlm.nih.gov/blast/db/ と blastn データベース ディレクトリ間の nt ファイル番号をクロスチェックします。
上記の bash スクリプトを使用して、不足している nt ファイルをダウンロードできます。
