解説
生体高分子を対象としたシミュレーション用構造の選択時の確認事項に関する指針
2024 年 23 巻 4 号 p. A20-A27
詳細
2024 年 23 巻 4 号 p. A20-A27
X線結晶構造などの実験的手法により構造が決定されProtein Data Bank(PDB)から利用可能になって以来，計算科学者はそれらを初期構造として，シミュレーション計算に用いるモデル構造を作成することは事実上必須になっていた．このモデル構造作成の過程は物理化学的に妥当なモデルを得るという分かりやすい目的のため，その手順は研究者によりさまざまであり，その作成方法はオープンソースプログラムを用いる方法から商用ソフトウエアを用いるなど多岐にわたる．一方で，その前処理のさらに前に複数の実験構造がPDBに登録されている場合，どの構造を優先して選ぶべきか，堅苦しくいうと「指針」のようなものは無いように思える．実際にFMOコンソーシアム等に参加する研究者や，シミュレーションを新規に試みる学生や研究者から相談を受ける事も多い．そこで，本稿ではシミュレーションを始める前に複数の構造を利用できる場合や1つの構造のみが利用できる時，どのような事に注目し分子ドッキング，分子動力学計算，フラグメント分子軌道法計算などのシミュレーションに用いるか，筆者の経験を元にして大まかに5つ程度の項目にまとめた．
Since the availability of experimentally determined structures, such as X-ray crystallography data from the Protein Data Bank (PDB), computational scientists have essentially been required to use these structures as initial models for simulations. The process of creating model structures varies widely among researchers, ranging from open source to commercial software, aimed at achieving physically and chemically accurate models. However, with multiple structures in the PDB, there is no clear guideline for selecting the optimal one. This article summarizes, based on the author's experience, five key points to consider when choosing structures for simulations like molecular docking, molecular dynamics, and fragment molecular orbital (FMO) calculations.
生体高分子を対象とする計算科学者は，座標データのダウンロードにProtein Data Bank (PDB)から公開されているタンパク質等のデータを使う機会は多いだろう．また，近年深層学習を用いた機械学習によるタンパク質の立体構造予測手法であるAlphaFold2 [1]によりほとんどの立体構造を予測することが可能になり，モデル構造を使用する機会はますます増えていくと思われる．実験構造もしくは予測モデルを高確率で利用可能となったため，それらの構造データを用いて座標データ以上の情報を得る，もしくは予測することが重要になってきた．タンパク質は主として20種のアミノ酸から構成されており，また非天然アミノ酸や糖鎖などの修飾をうけることがある．主に水素原子，炭素，窒素，酸素原子等から構成され，特に水素結合はαヘリックスやβシートなどの二次構造形成に寄与している．またトリプトファン，チロシン，フェニルアラニンは芳香環を含んでおり，CH/π相互作用やπ-π相互作用に関与する．これらの代表的な相互作用は好みのビューアーを用いて距離情報等を参考にして目視でも確認できるが，分子動力学法(MD)による動的構造解析やフラグメント分子軌道法(FMO)による電子的相互作用解析により豊富な情報を得る事ができる．
PDBに登録されている構造はそのままではシミュレーションに用いることはできない．よくある事例としては結晶構造には水素原子座標が存在しないため，適切な水素原子を発生させないとMD計算やFMO計算をかけることができない．またPDBの登録構造の中にも分子力場によるエネルギー値が不安定な部分構造が存在する場合もある．PDBからの座標公開は登録者の善意により成り立っている．よって計算科学者は自身のシミュレーションでうまく使えるような初期構造を自身の責任で妥当に前処理，つまりモデリングする必要がある．FMO法の提唱者である北浦先生はFu(Fu-site)というソフトウエアをご自身で作成している．(https://sourceforge.net/projects/fusuite/) 本ソフトウエアはGAMESSでFMO計算するための前処理をGUIから行うソフトウエアとなっている．前処理済みの構造は一般的なPDBファイルで保存可能であるので，その他のソフトウエアの入力ファイルとして使用することもできる．また理化学研究所のグループではAuto-FMO protocol [2]という手法を提案している．これはFMO計算の入力構造作成前処理ルーチンを自動化してFMO計算を行うプロトコールであり，内部でMOE [3]のSVLやPython script等を使って前処理全般を行う．また商用ソフトウエアに関しては国内海外含めて非常に数が多いため，筆者がある程度使った(大体10年くらい)経験を有するものに限ると，CCG社のMOE [3]のStructure PreparationモジュールやSchrödinger 社のMaestroから利用できるProtein Preparation Wizard [4]がある．どのソフトウエアも期待通りの動作を行ってくれるだろうから，これらソフトウエアを使用もしくは導入可能な環境にいるならばまずは試してみると良いだろう．
さらに，このようなモデリングによる構造準備のさらに前段階でよく質問される内容として，「複数の構造がある時はどれを選べばよいのか，またはどのあたりを注意すればよいか」がある．(Figure 1)こちらも大まかな基準はあるものの統一した見解はなくケースバイケースといえる．しかしこの際，本稿のスペースをお借りし，分子ドッキング，MD，FMOなどのシミュレーションに用いる構造選択の指針について，著者の経験に基づいた主観的項目（データ，目視，Viewer で確認すべき点）は以下のとおりである．まず確認する項目
If there are multiple structures to use in the simulation, we need to choose one.
1. 分解能が良いか(Sheldrick の基準(1.2Å) や 概ね2Å 前後の値)
2. 構造全般に欠損原子/ 残基がない，もしくは少ないか
3. 注目する部分構造の座標があるか
4. 注目部分に相互作用するリガンド分子があるか
5. 注目部分はクリスタルパッキングの影響を受けるか上の確認が終わったら状況に応じて確認する項目
a. R-free の値
b. 温度因子(B factor) の値
c. 注目している部分に電子密度が存在するか
d. Validation Report の各種項目
e. 水素原子の発生方法や形式電荷の検討
おそらく熟練の計算科学者の方はもっと他にもたくさんあると思うだろうが，本稿はこれからシミュレーションを始めようとする人に向けたものなので，わかりやすさ重視で5個までに抑えた．もっとある場合は各自で追加していただき，いらない場合は削除してもよい．また基本項目をおさえた後などの発展的な項目として，結晶構造の専門家と議論する際にしばしば話題に上がる指標なども後半に5つ加えた．本稿では主に前半の5つについて解説していく．
私はX線結晶構造の専門家ではないが，それでも分解能((Resolution)もしくは解像度)は知っており，構造ファイルを扱う上で最も有名な指標と言えるだろう．高い分解能(High resolution)は値が小さくなる．低い分解能(Low resolution)は値が大きくなる．例えば，分解能が1.0Åの構造は3.0Åの構造よりよい分解能であり，「低い分解能である」という表現は誤解を生じる．PDBのチュートリアルであるPDB-101のResolutionのセクションでも記述されているが「経験則として，高い分解能の構造は原子の位置がより信頼できる」 [5]という事になる．さて，ここで計算科学者がシミュレーションしたい1つのタンパク質種があって，そのタンパク質種について複数のPDB IDがPDBのWebサイトにあり，どれもX線結晶構造解析による構造であったとする．基本的に計算科学者は電子密度などの実験データから構造を構築する事はせず，提供されているPDBファイルやcifファイルをそのまま使いシミュレーション用のモデルを作るとする(余談だが私自身も最近結晶構造の専門家の方々から，PDBから電子密度データをダウンロードして自身でモデル構築する事を勧められるので，構築方法を勉強する必要があるかもしれない)．さて，この中からどれか1つを選ぶケースを考えると，まず分解能を比較する事が多いだろう．分解能についてはSheldrickの基準で「1.2Å」 [6, 7]というよく知られた数値があり，これは「原子分解能」と呼ばれることがあり，原子1つ1つを電子密度から区別することができるとされる．すべての実験構造の分解能が1.2Å以下ならばよいのだが実際はそういう事にはならない．
Figure 2は本稿執筆時のPDBにおける分解能の分布のヒストグラムであり，左に行くほど分解能の値が小さくなるため高分解能ということになる．この分布から1.8Åから2.0Åの範囲の構造が合計で31,361エントリと最も多く，1.2Å以下は5620エントリであり全20万エントリの2.7%程度に過ぎない．つまりSheldrickの基準をそのまま採用するならば，ほとんどの構造は，例えば2つの原子の位置を区別できないなどの，何かしらの「問題」をどこかに抱えていることになる．ところが実際にPDBからダウンロードしたファイルからシミュレーションに使った場合に問題が多かったかと自身で振り返ってみると，もちろんそうではなく著者は1.2Åに満たない場合でも十分に利用できた．1.2Åを閾値にするのは現実的ではないと考えられる．では分解能1.2Å以外ではどのような基準があるだろか．最近の報告ではDubachらのレビューでは一般に0.5~1Åの分解能はサブアトミック(または超高分解能)，1~1.5Å(Sheldrickの1.2Åを含む範囲)はアトミック，1.5~2Åは高分解能，2~3Åは中分解能，3~5Åは低分解能，5Å以上は超低分解能と述べている． [8]また既存のPDBのデータから再度モデルを構築する試みであるPDB REDOの開発に関する論文中で，入力ファイルの質の分類基準に1.2Å以下をアトミックと分類する [9]などのSheldrickの基準を意識しており，低分解能の基準についても記載がある．またDavid Blow著の「Outline of Crystallography for Biologists」 [10]の表11.1では分解能に応じて主鎖，側鎖，二次構造がどの程度判別できる大まかに記載されている．例えば，分解能が9Å程度でαヘリックス，4Å程度でβシート，2.4Å程度で芳香環の中心のくぼみ，1.5Å程度で原子が識別できる等である．この表の説明でも「分解能の解釈は注意して扱うべき」と述べられている [10]． また「Introduction to Protein Structure」 [11](注:日本語訳版あり)という著書では，分解能2 Åでロイシンとイソロイシン側鎖が区別できると記述されている．なお本稿の著者の私見ではだいたい2.0Å前後を目安にするものの，分解能の値による絞り込みはほどほどにして，そもそもターゲットの実験構造が存在するのか，さらに次の節で述べる自身の興味ある部分構造が存在するかを優先している．
Distribution of resolution of structures registered in the Protein Data Bank, data as of June 2024, total 207,040 entries(https://www.rcsb.org/stats/distribution-resolution)
ここまでに分解能に応じてどのような解釈をすべきか文献等を引用して紹介したがどうやら経験則の色合いが強い．分解能がよいに越したことはないが，それのみを絶対的基準とすることは避けるべきで，シミュレーションに用いる構造を自身でビューアーを用いて確認する事が特に重要だろう．また計算科学者がX線結晶構造解析を専門とする研究者と協力してシミュレーションを行う場合はより多くの指標について議題に上がる可能性が高い．これはリストの下の行に加えた．例えば分解能以外の指標としてR-free [12]や温度因子(B factor) [13]があり，PDBのWebサイトから構造の評価であるvalidation report [14]が各PDB IDについて利用できる．また電子密度を描画することで原子の位置が信頼できるか確認してから判断する事もある．より多くの情報から判断する必要がある読者はそれらを参考にしてほしい．
ここで残りの指針を解説する前に，図3に著者がよくみかける構造のパターンについて5つほどまとめた．この図では円をタンパク質などの高分子とし三角形の切れ込みが興味のある結合サイトを表す受容体分子(Receptor, Rと表記)とリガンド分子(Ligand, Lと表記)の相互作用を表しており，この分野で受容体分子とリガンド分子の相互作用の表現やいわゆる"lock-and-key model"を説明する時しばしば使われる伝統的なモデルである．研究者によってはここで示した5つ以外の状態を定義しているが，本稿はシミュレーション初心者に向けているのでここでは触れない．図3.Aは受容体分子の結合サイトに対して三角形のリガンド分子が結合している状態を示している．なお受容体はタンパク質，リガンド分子は低分子を想定するが，この概念はペプチド分子などの中分子，タンパク質間相互作用における高分子にも拡張できるため相対的に決定される．もし自身のシミュレーションでタンパク質とリガンド分子の相互作用解析を行う場合，原子欠損等がなく高分解能で図3.Aのような状態のリガンド分子結合構造を入手できる場合は，まずは優先して選択するのがよい．Figure 3.BからEについての詳細は次節以降で解説する．
Example of receptor (R) and ligand (L) interaction using a model with a circular triangular cutout (A) General R-L complex structure (B) There is an invisible substructure at the binding site (C) L is not present (D) The structure of the binding site changes significantly in the absence of L (E) The binding site changes due to crystal packing
実験構造については，座標のどこかに原子や残基の欠損がある場合がある．これは実際にその部分にアミノ酸残基などが存在しない場合もあるだろうが，分子中に残基は存在するが結晶構造データから座標が構築できないケース(いわゆるmissing residueやatoms)がある．ここでは後者について取り扱う．ところで原子や残基に欠損があるとシミュレーションに何故問題になるのだろうか?一般的には，シミュレーションは何かしらのエネルギーやスコアを計算し，エネルギー関数は(カットオフを導入しない場合)系の全原子について計算するため，原子が存在しないとその分だけ不正確になる．量子化学計算では正常にシミュレーションできないばかりではなく，計算が収束したとしても想定している分子と異なるため意図したシミュレーションができないことになる．MD計算では力場パラメータを適切に割り振れない可能性が非常に高い．よって原子の数を過不足なく整える必要があるので欠損部分をホモロジーモデリングや適切なループの移植で補うか，補完が難しいなら適切に末端をキャップし(アセチル基，ACEと表記もしくはN-メチル基，NMEと表記)意図しない静電相互作用が発生しないようにする必要がある．現時点では計算結果の使用用途に影響を及ぼさない簡単に使えるフリーな構造補完手法はあまりないが，例えば最近だと簡便に利用できるものとしてはWebサーバー,Protein Repair and Analysis Serverのソースコードが存在する [15]．
シミュレーションをする際にはタンパク質の全体のフォールドの動的挙動の解析もあるが，興味ある特定の部分構造に注目するケースが多いだろう．例えば図3.Aで示したような創薬のターゲットタンパク質の結合サイトやタンパク質間の相互作用部位など，広い意味で受容体分子とリガンド分子の相互作用がある．注目している結合サイトの構造が存在し，遠い部分にのみ欠損残基や原子がある場合は，用いるシミュレーションの種類により状況が異なる．例えば分子ドッキングにより計算されるドッキングスコアは，MDやFMOに比べるとスコア計算は簡略化されているので，得られた予測結合状態が相互作用解析の役に立つ事が期待できるが，MD計算やFMO計算では前述のように少量の欠損も許容されない事がほとんどだろう．ただし，ドッキング計算においても結合サイトの近くに原子欠損がある場合はスコアやエネルギーが正確に計算できない可能性が高いので結果をそのまま利用する事は避けた方がよい．そのような構造しか選ぶことができない状況ならば，適切な構造補完手法で座標を用意しなければならない．
受容体の原子がすべてあるからといって，そのままのすべての構造をシミュレーションに含めてはいけない場合もある．例えばPDBの登録構造には結晶化のために人工的にアミノ酸配列を置換している場合等があり，おそらく変異前の構造のシミュレーションに興味がある事が多いので，元の残基に戻す方が適切である．このような変異残基に関する情報は原著論文やPDBファイルに記載がある事が多い．一方で，翻訳後修飾等によるリン酸化されたアミノ酸残基や糖鎖修飾構造など，生体中構造としては正しいのだが力場パラメータを設定できずMD計算等を実行できない場合もある．このような時はシミュレーション自体をあきらめるのではなく，そのような残基を力場パラメータが存在する天然アミノ酸等に戻す，または削除した上でシミュレーションし，まずは天然アミノ酸での挙動を研究する事も選択肢として考えられる．
ターゲットタンパク質の興味のある部分に座標が存在した場合，次に注目すべきは相互作用相手となる分子の座標が結合サイトに存在するかである．Figure 3.Cのような受容体中にリガンド分子の座標が存在しない場合がある．このような構造しか入手できない場合はそのままの状態を用いてシミュレーションを行う事になるが，一般的には何かしら相互作用する分子が存在する事が好ましい．なぜならリガンド分子が存在しない事で結合サイトがいわば「閉じた」状態になり，相互作用に適した形ではなくなっている可能性が考えられるからだ．例えばAlphaFold2の構造に対するリガンドドッキングのベンチマークに関する論文 [16]ではリガンド分子が結合していない状態では，ドッキングプログラムによる予測性能が減少する傾向があることが示されている．これは結合サイトの主鎖や側鎖コンフォメーションがリガンド分子との相互作用に適していない場合，たとえドッキングプログラムが計算過程で正解のコンフォメーションを予測できていたとしても，結合サイトが閉じているならば原子間衝突がおこりスコアが低くなるためと考えられる．このような「閉じた」構造から「開いた」構造を得る目的ではMD計算によるコンフォメーションのサンプリングが有力である．もし複数の実験構造を利用できる場合は結合サイトにリガンド分子が結合しているかを確認する事が重要となる．
なにかしらの低分子があるからといって，そのまますべての構造をシミュレーションに含めてはいけない場合もある．登録構造には生体内に存在しない添加剤が含まれる事がある．(有名なのはPEGやSO4，生体機能に関係ない金属イオンなど) むしろシミュレーションが終わった後によくない構造を使っていた事がわかると，最悪もう一回計算をやり直すことになるため，もし結晶構造の研究者と共同で研究している場合はそれらの分子をシミュレーションに含めるべきか相談する事を勧める．
以上で述べた項目はPDBファイルをビューアーで読み込み目視により直接確認することができる．一方で結晶構造の場合，PDBファイル中に見かけ上分子が1つしか存在しなくても，結晶中では規則正しく整列して存在し結晶中では隣接する分子が存在する．これはクリスタルパッキング(crystal packing)と呼ばれる [17]．結晶構造中で隣り合っている分子が，生体内中および水溶液中でそのままの相対配置を保っているとは限らない．さらに気を付けるべきは図3.Eのようにクリスタルパッキングにより注目している結合サイトの構造が変化してしまっている可能性があることだ．PDBファイルやcifファイルではその構造が生物学的に機能する最小単位なのか，クリスタルパッキングなのか調べるための情報を提供している．17PDB形式の場合，BIOMTやMTRIXというレコードに存在するが，各種ビューアーにはこれらの情報から結晶中の複合体構造を生成する機能をもつためそれらを利用するとよい．例えば著者らのグループがよく使っているMOEではPDBファイルを読み込むときに生物学的複合体を生成する選択肢を選べることができる．(図4) Symmetryの項目からBiomolecule Assemblyを選ぶと生物学的単位を生成できる．発生する構造は複合体もしくは単体のどちらの場合もある．またContactsが含まれている項目を選択すると結晶中での相互作用を発生できるので，クリスタルパッキングの可能性の検討を目視で確認できる．例に示したH-PGDSと阻害剤の複合体構造 [18]で本操作を行ったところ，結合サイト周辺には他の結晶中構造が少なくクリスタルパッキング影響は少ないと判断できる．もしクリスタルパッキングの可能性が高い場合は，そのようなリスクを考慮した上でそのまま用いるか，MD計算などで構造を緩和してから使用する事を検討する必要がある．
An example of generating a crystal structure from the window displayed when a PDB file is loaded in MOE. In this example, after downloading PDB ID: 5YWX, "Biomolecule+Contacts" is used to display the residues (green) with atoms at 4.5 Å from the biological unit. The CPK model in the center is the ligand molecule.
本稿では，生体高分子を対象とした分子ドッキング，MD，FMOなどのシミュレーションをこれから試みようとする学生や研究者にむけて，シミュレーション前に行う構造選択の指針を5つの点に絞り，その理由や具体的な例を交えて説明した．まず座標ファイルの全体的な質を調べるためには分解能を確認する事が重要であるが，分解能の絶対的な閾値はない．使用するシミュレーションの種類や，ターゲットタンパク質の全体に注目するのか，部分構造に注目するかで対応が変わる．結合サイトを確認する際には相手側分子の有無により結合様式に違いが出る可能性が高いので，何かしら分子が結合している方が望ましい．結晶構造を目視する上ではPDBファイルには座標が存在しなくても結晶中には近くに分子が存在する事がありクリスタルパッキングによる影響の有無を確認するのがよい．
近年PDBに登録される構造はクライオ電顕構造を含めて増加傾向にあり，今後は計算機によるシミュレーションはますます利用されていくことが想定されるが，クライオ電顕による構造も本稿で参照した基準は概ね当てはまると期待できる．一方でNMRにより得られる構造は1つのPDB IDに複数のモデルが登録されている事が多く各モデルを注意深く選ぶ必要があるが，本稿の注目すべきポイントはある程度参考になるだろう．
実験構造からシミュレーションに用いる構造を準備する場合は，ここで述べた優先的な5項目に注目して構造を目視で確認しながら，妥当なシミュレーション条件を検討し，状況に応じてR-freeやB factor，電子密度を確認してほしい．そして，シミュレーションの結果の解析段階，成果発表の直前，さらには発表中に入力構造の不備を他の研究者から指摘されないように努めてほしい．
本特集号の執筆に招待していいただいた大阪大学大学院薬学研究科量子生命情報薬学分野の福澤薫先生に深く感謝いたします．
