mtDNAの系統解析などをする際に、同じ塩基配列(ハプロタイプ)をもつサンプルを取り除きたいことがある。Fastaファイルから重複した配列を手作業で取り除くのは手間がかかるので、seqkitを使って一瞬で終わらせよう。
追記:Twitterでvsearchが便利ということを教えていただきました。同一配列の除外と検出については、完全にseqkitの上位互換ですね。情報提供ありがとうございます!
環境
macOS 12.3.1
seqkit v2.0.0
vsearch 1.2.11
- Macの場合はHomebrewで導入できる
brew install seqkit vsearch
- condaは
conda install -c bioconda seqkit
同一の配列を除去する
seqkitでやる場合。
ターミナル
seqkit rmdup -s samples.fasta > out.fasta
# 同じ名前の配列を除去する場合は-sではなく-nを使う
[INFO] 245 duplicated records removed
注意
同じハプロタイプでも、塩基配列長が違うと除去されないので注意。
vsearchでやる場合。
ターミナル
vsearch --cluster_fast samples.fasta -id 1 --centroids out.fasta
Reading file samples.fasta 100%
300049 nt in 472 seqs, min 526, max 655, avg 636
Masking 100%
Sorting by length 100%
Counting k-mers 100%
Clustering 100%
Sorting clusters 100%
Writing clusters 100%
Clusters: 153 Size min 1, max 32, avg 3.1
Singletons: 94, 19.9% of seqs, 61.4% of clusters
注意
-のようなGAPが配列中にあるとエラーになるので、予め除去しておく(アライメント不要)。
塩基配列長が違っても一致した配列は除去される。
簡単!
同一の配列をもつサンプル名のリストやハプロタイプを取得する
同一配列の除外ができるということは、ハプロタイプの検出にも使える。
ここではvsearchを使う。
#ハプロタイプごとに分割されたFastaファイルの出力
vsearch --cluster_fast samples.fasta -id 1 --clusters hap
#各サンプルとハプロタイプの対応リストの出力
vsearch --cluster_fast samples.fasta -id 1 --mothur_shared_out list.tsv
便利!
注意
サンプル名にハイフン(-)が入っていると、その後の文字列がうまく認識されないので、アンダーバー(_)などに置換しておきましょう。
参考サイト