HMMERとは
HMMERは、タンパク質やDNA配列のドメインファミリーを確率論的にモデル化するツールを提供するソフトウェアパッケージです。この記事では、HMMERで使用できるプログラムの使い方を紹介します。
動作環境
- Windows11
- Ubuntu 22.04.3 LTS
インストール方法
HMMERはanacondaでインストールできます。(他にも方法はあります)
anacondaの環境構築方法は他の記事でまとめています。
$ conda install bioconda::hmmer
HMMER プログラム一覧(抜粋)
HMMERで使用できるプログラムを抜粋しました。他のプログラムは、ユーザーガイドから使用方法を確認できます。
1. hmmstat
2. hmmpress
3. hmmscan
4. hmmfetch
5. hmmsearch
これらプログラムの使用方法の紹介にあたり、HMM profiles, 2種のアミノ酸配列, 大腸菌のゲノム配列を用意しました(以下データセット)。
データセット
HMM profiles
HMM profilesはInterPro(Pfam)にあるものをダウンロードしました(Pfam-A.hmm.gz)。ファイルは圧縮されているので$ gunzip Pfam-A.hmm.gz
で解凍する必要があります(Pfam-A.hmm)。
アミノ酸配列
アミノ酸配列を二つ用意しました。これらのアミノ酸は、ドメインを一つまたは複数持っています。
>singledomain
MSKSDVFHLGLTKNDLQGATLAIVPGDPDRVEKIAALMDKPVKLASHREFTTWRAELDGKPVIVCSTGIGGPSTSIAVEELAQLGIRTFLRIGTTGAIQPHINVGDVLVTTASVRLDGASLHFAPLEFPAVADFECTTALVEAAKSIGATTHVGVTASSDTFYPGQERYDTYSGRVVRHFKGSMEEWQAMGVMNYEMESATLLTMCASQGLRAGMVAGVIVNRTQQEIPNAETMKQTESHAVKIVVEAARRLL
>multipledomain
MADSPLRHLSTFVRYLRGERTAVALLALLVPAGMGLQLVAPYLLRGFIDGALSGDSRKTLLDLAAWSLAAAVGTLVVTAGTEALSSRVAWRSTNRLRADLVEHCLSRPPGFYRKHPPGELVERMDGDVTRLAAVMSTLLLELLAQALLIVGILVALFRLEWRLALVFAPFAAGTLLLLRTLVGRAMPFVTARQRVAADLQGFLEERLAAAEDLRVNGASRYTLRELGDRQDDLYRKARDAARASVRWPATVQGLSAVSVVLALAVSAWLHARGQLSTGTAFASLSYAMLLRRPLLAVTTRFRELEDAAASAQRLRDLLGHGTAAPRTGRGTLPAGLPGVRFDGVSFGYEPDEPVLRDVSFTLRPGERLGVVGRTGSGKSTVVRLLFGLHHPGAGSVSAGGLDLTEIDPRALRSRVALVTQEVHVFYASLRDNLTFFDRSVPDDRLRAALGEAGLGPWLRTLPDGLDTPLGAGARGMSAGEEQQLALARVFLRDPGLVLMDEPTARLDPYSERLLMPALERLLEGRTAVVVEHRPHLLRNVDRILVLEEGKVAEEGERRVLAADPGSRFHALLRTAGATR
ゲノム配列
ゲノム配列は、NCBIから取得しました。今回使用するのは、Escherichia coli str. K-12 substr. MG1655の配列です。以下スクショのように指定してダウンロードしました。
GCA_000005845.2のフォルダにあるprotein.faaというファイル名を、GCA_000005845.2.faaに変更して使用しました。
実施例
1. hmmstat
hmmstatは、今回使用するPfam-A.hmm(hmmfile)の要約統計量を出力します。これにより、HMM profilesに含まれているプロファイルが簡単に分かります。
$ hmmstat Pfam-A.hmm > pfam_hmmstat.txt
$ head -n 15 pfam_hmmstat.txt
# hmmstat :: display summary statistics for a profile file
# HMMER 3.1b2 (February 2015); http://hmmer.org/
# Copyright (C) 2015 Howard Hughes Medical Institute.
# Freely distributed under the GNU General Public License (GPLv3).
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
# - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
#
# idx name accession nseq eff_nseq M relent info p relE compKL
# ---- -------------------- ------------ -------- -------- ------ ------ ------ ------ ------
1 1-cysPrx_C PF10417.14 30 15.00 41 1.34 1.28 1.26 0.12
2 10_blade PF21734.2 14 1.62 345 0.59 0.55 0.54 0.07
3 117-like_vir PF21578.2 1 0.38 101 0.59 0.54 0.51 0.08
4 120_Rick_ant PF12574.13 2 0.45 240 0.59 0.60 0.52 0.03
5 12TM_1 PF09847.14 6 2.00 448 0.59 0.62 0.51 0.12
2. hmmpress
hmmpressは標準的なhmmfileをバイナリ圧縮します。これは、hmmscanを実施する前に必要なステップです。出力として、hmmfile.h3m, hmmfile.h3m, hmmfile.h3f, hmmfile.h3pが得られます。実際にこのコマンドを使用すると以下の出力が得られます。
$ hmmpress Pfam-A.hmm
Working... done.
Pressed and indexed 21979 HMMs (21979 names and 21979 accessions).
Models pressed into binary file: Pfam-A.hmm.h3m
SSI index for binary model file: Pfam-A.hmm.h3i
Profiles (MSV part) pressed into: Pfam-A.hmm.h3f
Profiles (remainder) pressed into: Pfam-A.hmm.h3p
3. hmmscan
hmmscanは、クエリー配列のファイル(今回はsample.fa)にある各配列をもとに、hmmdb(hmmpressしたhmmfile)を検索し、配列に最もマッチしたプロファイルのリストを出力します。以下出力結果から、単一ドメインの配列には1つのプロファイルがヒットし、それはPNP_UDP_1であることがわかりました。複数ドメインの配列には複数のプロファイルがヒットしました。
$ hmmscan --domtblout sample_hmmscan.txt Pfam-A.hmm sample.fa
$ head sample_hmmscan.txt #出力確認
# --- full sequence --- -------------- this domain ------------- hmm coord ali coord env coord
# target name accession tlen query name accession qlen E-value score bias # of c-Evalue i-Evalue score bias from to from to from to acc description of target
#------------------- ---------- ----- -------------------- ---------- ----- --------- ------ ----- --- --- --------- --------- ------ ----- ----- ----- ----- ----- ----- ----- ---- ---------------------
PNP_UDP_1 PF01048.25 233 singledomain - 253 4e-40 138.0 1.5 1 1 2e-44 4.4e-40 137.8 1.5 1 223 20 244 20 253 0.85 Phosphorylase superfamily
ABC_membrane PF00664.28 274 multipledomain - 579 1.1e-29 104.3 0.0 1 1 5.3e-33 1.3e-29 104.0 0.0 13 274 35 294 23 294 0.96 ABC transporter transmembrane region
ABC_tran PF00005.32 137 multipledomain - 579 3.4e-27 96.0 0.0 1 1 2.3e-30 5.6e-27 95.3 0.0 1 137 355 504 355 504 0.95 ABC transporter
AAA_21 PF13304.11 304 multipledomain - 579 5.1e-06 27.0 0.5 1 2 0.0037 8.9 6.5 0.0 3 22 369 388 368 414 0.78 AAA domain, putative AbiEii toxin, Type IV TA system
AAA_21 PF13304.11 304 multipledomain - 579 5.1e-06 27.0 0.5 2 2 7.4e-07 0.0018 18.7 0.1 236 295 474 530 421 537 0.88 AAA domain, putative AbiEii toxin, Type IV TA system
DUF87 PF01935.22 220 multipledomain - 579 0.0074 16.7 0.1 1 1 6.4e-06 0.016 15.6 0.1 25 44 367 386 361 388 0.92 Helicase HerA, central domain
SMC_N PF02463.24 220 multipledomain - 579 0.012 15.3 0.0 1 2 0.01 26 4.4 0.0 28 43 369 384 357 395 0.84 RecF/RecN/SMC N terminal domain
4. hmmfetch
hmmfetchは、複数のHMMプロファイルがまとまったhmmfileから、目的とするHMMプロファイルを抽出します。ここでは、単一ドメインの配列にヒットした、PNP_UDP_1のプロファイルを取り出します。
$ hmmfetch -o pnp_udp_1.hmm Pfam-A.hmm PNP_UDP_1
$ head -n 5 pnp_udp_1.hmm #出力確認
HMMER3/f [3.1b2 | February 2015]
NAME PNP_UDP_1
ACC PF01048.25
DESC Phosphorylase superfamily
LENG 233
$ tail -n 5 pnp_udp_1.hmm #出力確認
0.02583 6.16670 3.75488 0.61958 0.77255 0.48576 0.95510
233 1.95585 6.02370 2.96175 2.00713 5.37000 3.63762 3.97967 3.70400 2.10203 3.54825 5.06404 2.88035 3.47748 3.13431 2.06920 2.33703 2.85980 3.92648 6.45565 5.04395 491 a - - H
2.68618 4.42225 2.77519 2.73123 3.46354 2.40513 3.72494 3.29354 2.67741 2.69355 4.24690 2.90347 2.73739 3.18146 2.89801 2.37887 2.77519 2.98518 4.58477 3.61503
0.00215 6.14302 * 0.61958 0.77255 0.00000 *
//
5. hmmsearch
hmmsearchは、配列データベースに対して、目的のプロファイルにマッチする配列があるかを検索することができます。hmmfetchで取り出したPNP_UDP_1のプロファイルにマッチする配列が、Escherichia coli str. K-12 substr. MG1655に存在するか調べました。その結果、複数の配列にこのドメインが含まれていることがわかりました。
$ hmmsearch --domtblout GCA_000005845.2.txt pnp_udp_1.hmm GCA_000005845.2.faa
$ head -n 8 GCA_000005845.2.txt
# --- full sequence --- -------------- this domain ------------- hmm coord ali coord env coord
# target name accession tlen query name accession qlen E-value score bias # of c-Evalue i-Evalue score bias from to from to from to acc description of target
#------------------- ---------- ----- -------------------- ---------- ----- --------- ------ ----- --- --- --------- --------- ------ ----- ----- ----- ----- ----- ----- ----- ---- ---------------------
AAC73270.1 - 232 PNP_UDP_1 PF01048.25 233 8.2e-57 190.2 2.2 1 1 1.3e-59 9.1e-57 190.1 2.2 1 232 2 226 2 227 0.96 5'-methylthioadenosine/S-adenosylhomocysteine nucleosidase [Escherichia coli str. K-12 substr. MG1655]
AAC75460.1 - 277 PNP_UDP_1 PF01048.25 233 7e-45 151.2 0.6 1 1 1.1e-47 7.9e-45 151.0 0.6 1 233 27 274 27 274 0.95 xanthosine phosphorylase [Escherichia coli str. K-12 substr. MG1655]
AAC77337.1 - 239 PNP_UDP_1 PF01048.25 233 4.8e-41 138.6 1.3 1 1 7.8e-44 5.6e-41 138.4 1.3 2 218 16 223 15 234 0.90 purine nucleoside phosphorylase [Escherichia coli str. K-12 substr. MG1655]
AAC76834.1 - 253 PNP_UDP_1 PF01048.25 233 7.7e-41 138.0 1.5 1 1 1.2e-43 8.7e-41 137.8 1.5 1 223 20 244 20 253 0.85 uridine phosphorylase [Escherichia coli str. K-12 substr. MG1655]
AAC75046.1 - 484 PNP_UDP_1 PF01048.25 233 3.5e-22 76.9 0.0 1 1 8.1e-25 5.8e-22 76.2 0.0 37 201 266 431 246 452 0.86 AMP nucleosidase [Escherichia coli str. K-12 substr. MG1655]
まとめ
HMMERの使用方法を簡単にまとめました。以降の記事で、詳細なオプションの紹介や、出力ファイルの説明等していきたいです。