JP4309282B2 - 複数鎖を有するタンパク質の立体構造構築方法 - Google Patents

Description

技術分野
本発明は複数鎖を有するタンパク質の立体構造構築方法、詳しくは立体構造が未知で、複数鎖のタンパク質の立体構造を予測する方法に関する。この方法によれば、複数鎖より成るタンパク質を単純化のために擬似単数鎖のように扱い、当該タンパク質を構成する複数のペプチド鎖間の相互作用を考慮した構造予測を行うことができる。これにより、後述の実施例で述べるように、従来法に比較して、より信頼性の高いタンパク質立体構造予測を行うことができる。
更に、本発明は、この方法で得られるタンパク質の立体構造モデル、前記方法に使用可能なデータベース、データベース構造、及びコンピュータソフトプログラム、これを搭載したコンピュータやインターフェース等にも関する。
例えば、複数鎖から成るアミノ酸配列の各Ｃ末端残基に区切り文字である終止記号を付加し、主鎖と側鎖の３次元座標から成るタンパク質の立体構造を構築する際、この終止記号から各鎖Ｃ末端の残基番号を得ることで、タンパク質鎖間のアミノ酸残基間の相互作用を考慮したポテンシャル計算を行うことができる。この方法により、側鎖のパッキングが優れているモデルが構築できるので、従来法に比較して、信頼性の高い複数鎖より成るタンパク質立体構造予測を行うことができる。尚、区切り文字の付加は、複数の鎖を全て相互に結合して単数鎖として扱えるようにするものであり、単数鎖にしたときに結合する最後の鎖のＣ末端部（単数鎖にしたときの最末端部）には必要がない。
本発明において、ポリペプチド鎖、リガンド等、タンパク質を構成する、又は構成することのできる鎖、即ち当該タンパク質の鎖を、単に「タンパク質鎖」と称することがある。
タンパク質中にポリペプチド鎖以外の鎖（タンパク質鎖）が存在し、例えば複数鎖を構成する一つのタンパク質鎖であるコンポーネント（成分）がペプチド等の低分子リガンドである場合には、そのリガンドのアミノ酸配列を任意に改変することで多様な結合状態を作成することができる。
また、ポテンシャルパラメータとして、アミノ酸残基を構成する原子間のパラメータに、任意の原子間のポテンシャルパラメータを追加することで、任意のリガンド分子への改変を行うことができる。更に、リガンドのアミノ酸配列を固定し、周囲のタンパク質鎖のアミノ酸配列を多様に改変したデータセットを作成することで、特定のリガンドに結合し得る多様な受容体モデルを構築することができる。また、複数鎖を構成しているタンパク質間の相互作用を正確に表現できるので、機能と関係した認識部位を記述したモデルを構築することができ、相互作用部位のアミノ酸残基を改変することによって、機能の増減を調節できるようなモデルを構築することができる。更に、単数鎖の場合において、領域（ドメイン）や部品（モジュール）に切り分けることにより複数鎖と見なした後、これを擬似単数鎖として復元することで立体構造モデルの高精度化を図ることができる。
尚、本発明のタンパク質の立体構造構築方法において、その基本骨格は、タンパク質の経験的なモデリング法、特にホモロジーモデリング法或いはスレッディング法を利用するものであり、その際Ｘ線結晶構造解析等により、タンパク質の立体構造として３次元座標が決定されたものを参照すべき鋳型に用いる。このとき、実質的には、複数鎖として構造決定されたタンパク質を鋳型として用いなければ、各タンパク質鎖間、特に各ポリペプチド鎖間の正確な相対配置は決められないことが多い。換言すれば、本発明は、相対配置が既知である参照タンパク質複合体を用いる経験的なモデリング法である。但し、モデリングの対象となる目的タンパク質と参照タンパク質の間で、アミノ酸配列のホモロジーは、必ずしも高いことを必要とせず、後述する規定の条件（Ｅ−ｖａｌｕｅ）を満たす全てのタンパク質を参照タンパク質として用いることができる（スレッディング（ｔｈｒｅａｄｉｎｇ）法）。また、例えば相互作用界面に関してのみモデリングするので良ければ、必ずしも複合体全体について参照し得る実験構造が存在することを要しない。
背景技術
世界中でゲノムの配列の決定が行われてきており、７０種を超えるゲノムのコードするアミノ酸配列も同様に多数決定されてきている。特に、ヒトを頂点とし、多くの生物種において、全ゲノム配列解析が進められて、その配列情報がデータベース化されている（Ｇｅｒａｒｄｏ Ｊｉｍｅｎｅｚ−Ｓａｎｃｈｅｚ，Ｎａｔｕｒｅ ４０９，８５３−８５５（２００１）参照。）。ゲノム配列から遺伝子の機能を特定、予測することがある程度可能であるが、配列情報のみでは機能を予測できない遺伝子が非常に多く存在する。遺伝子は、そのＤＮＡ配列からアミノ酸配列に翻訳されたタンパク質が実際に機能を担っている。タンパク質の機能解明には、Ｘ線結晶構造解析、ＮＭＲ等の実験的手法により、その立体構造を決定することが行われているが、一般に相当量の時間と労力を要する作業を必要とする。そこで、目的とするタンパク質とアミノ酸配列の相同性の高い立体構造既知のタンパク質を探し出し、目的タンパク質の機能を予測することが行われている。
立体構造未知のタンパク質アミノ酸配列から、単体としてのタンパク質立体構造を予測する方法としては、ホモロジーモデリング法が一般的である（Ｔ．Ｙｏｎｅｄａ，Ｈ．Ｋｏｍｏｏｋａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊ．Ｐｒｏｔｅｉｎ Ｃｈｅｍ．，１６，５９７−６０５，１９９７参照。この全内容は参考により本明細書中にその一部として組み込まれる。）。これは、主に次の４つのステップより成る計算科学的手法である。
▲１▼ 任意の立体構造が未知の目的タンパク質のアミノ酸配列（目的配列）が与えられたとき、目的配列と類似の配列を有する参照タンパク質（参照配列）をＰＤＢ（プロテインデータバンク）のような立体構造データベースから、配列が統計的に有意に似ているものを検索（ホモロジー検索）し、目的配列と参照配列の間で配列を並置したもの（アライメント）を与える。
このデータベース検索とアライメントを行うためには、ＦＡＳＴＡ（Ｐｅａｒｓｏｎ ＷＲ，Ｍｅｔｈｏｄｓ Ｅｎｚｙｍｏｌ，２６６，２２７−２５８，１９９６参照。）、ＰＳＩ−ＢＬＡＳＴ（Ｓｃｈａｆｆｅｒ ＡＡ，Ｗｏｌｆ ＹＩ，Ｐｏｎｔｉｎｇ ＣＰ，Ｋｏｏｎｉｎ ＥＶ，Ａｒａｖｉｎｄ Ｌ ａｎｄ Ａｌｔｓｃｈｕｌ ＳＦ，Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ，１２，１０００−１０１１，１９９９参照。）、ＬＩＢＲＡ（Ｏｔａ，Ｍ．ａｎｄ Ｎｉｓｈｉｋａｗａ，Ｋ．，Ｐｒｏｔｅｉｎ Ｅｎｇｉｎｅｅｒｉｎｇ，１０，３３９−３５１，１９９７参照。）等のコンピュータソフトを使用することができる。ＦＡＳＴＡは２０種類の天然アミノ酸を意味する２０種のアルファベット文字配列のマッチングを行うプログラムであり、高ホモロジー（アミノ酸残基の一致度約３０％以上、ＦＡＳＴＡのＥ−ｖａｌｕｅでは約０．０１以下に相当。）の参照タンパク質に対して立体構造構築をすると、信頼性の高いモデルを構築できると報告されている。
一方、ＰＳＩ−ＢＬＡＳＴでは、同じように文字配列のマッチングを行うが、文字が一致しているか否かの情報ではなく、プロファイルと呼ばれる文字の一致の度合いを類縁タンパク質の文字配列部位ごとの置換行列として算出し、更に繰り返し計算を行うことによりアライメントを最適化する性質を持っている。また、ＬＩＢＲＡは３Ｄ−１Ｄ法（スレッディング法の別名）に基づくプログラムであり、既知立体構造を精査し、目的配列に対して類似配列を検索するため、ＦＡＳＴＡやＰＳＩ−ＢＬＡＳＴとは検索アルゴリズムが明らかに異なる。そのため、ＦＡＳＴＡやＰＳＩ−ＢＬＡＳＴとは異なって間違いも含むが幅広く配列間類似性を指摘できる場合がある。
▲２▼ ＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ、ＬＩＢＲＡ等により算出したアライメントを用いれば、目的配列と類似した参照配列間の文字配列部位ごとの対応関係が決まる。この関係に基づき、参照タンパク質の３次元座標から目的配列のアミノ酸残基の３次元座標を作成する。
▲３▼ 参照配列に対し目的配列側に対応するアミノ酸残基が存在しない場合（アミノ酸残基欠損）には、参照タンパク質側のその位置のアミノ酸残基座標は用いず、逆に、目的配列に対し参照配列側に対応するアミノ酸残基が存在しない場合（アミノ酸残基挿入）には、その位置の目的配列のアミノ酸残基座標は、予め用意しておいたタンパク質断片座標データベース（例えば、Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。この全内容は参考により本明細書中にその一部として組み込まれる。）から適切なものを検索して作成する。
▲４▼ 上記▲２▼及び▲３▼によるタンパク質座標の構築では、アミノ酸残基間に構造的に不適切な隙間や衝突や歪みが生じることがあるので、エネルギー極小化計算や分子動力学計算により、これらの構造的な歪みを解消する（Ｍ．Ｔａｋｅｄａ−Ｓｈｉｔａｋａ，Ｈ．Ｕｍｅｙａｍａ，ＦＥＢＳ Ｌｅｔｔｅｒｓ，４２５，４４８−４５２，１９９８参照。この全内容は参考により本明細書中にその一部として組み込まれる。）。
モデリングソフトによっては、▲４▼の構造的な歪みの解消をスムーズに行うため、上記▲２▼〜▲４▼の計算及び検索処理をタンパク質全原子に対して同時に行うのではなく、例えばシミュレーティッドアニーリング（ＳＡ）法等を用いて、段階的に行うものもある。
「ＳＡ法」とは、或る系の状態ｘに対して摂動を加え、新しい状態ｘ’を得たとき、新しい状態でのエネルギー値Ｅ（ｘ’）が旧状態のエネルギー値Ｅ（ｘ）より小さければ高い確率で、大きければ低い確率で新状態ｘ’への遷移を行うことにより、局所最小点に捉えられることなくエネルギーＥの大域的最小点を発見しようというものである。即ち、先ずタンパク質の骨格を形成するα炭素原子について、次いでα炭素原子を含む主鎖原子について、ＳＡ法によるエネルギー極小化を行い、最後に側鎖原子を含むタンパク質原子全体についてＳＡ法によるエネルギー極小化を行うものである。「分子力学計算」とは、系のポテンシャルエネルギーＥを座標の関数として表し、主に最大傾斜法又は共役勾配法等のエネルギー極小化計算によって、ポテンシャルエネルギーＥの下がるべき座標変化を計算する方法である。「モンテカルロ法」とは、統計力学に基いた確率論的エネルギー最適化計算法である。
以上により、▲２▼の目的配列に対するアライメントが得られれば、３次元座標の作成を通じてその立体構造を予測構築することができる（前記Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。
しかしながら、従来は、複数鎖を有するタンパク質に対する任意の配列が与えられたとき、ＰＤＢのような立体構造データベースから、各鎖に対して独立に、類似性の高いタンパク質を単数若しくは複数選び出し、アライメントを与え、それぞれ独立にモデリングするため、後述の実施例において説明するように、タンパク質−タンパク質間の相互作用様式を十分に反映した立体構造を得ることができなかった。特に、医薬品や酵素の基質等のリガンド結合部位が複数鎖から構成される場合には、上記の欠点は重大である。また、反対に医薬品や酵素の基質に対して、それ等リガンドを受容するタンパク質の機能効率等の性質を改変する場合に、同様の欠点が発生する。従来、この欠点を取り除く方法としては、分子動力学法等により、モデル構築後に、複数鎖間の相互作用様式を推定することで改良する方法があるが、グローバルミニマムを探索するためには膨大な計算時間と計算機資源が必要なため、また、本来多体に対応すべき分子力学パラメータが２体対応のため真の立体構造を与え難いために、ゲノムワイドな産業利用には適さない等の問題点がある。
発明の開示
１．発明が解決しようとする課題
本発明者等は、複数鎖を有するタンパク質に対する任意の配列が与えられたとき、当該複数鎖を単数鎖に変換してその立体構造を構築する方法について検討を行った。実験的には、天然又は概念上の複数鎖を単数鎖に変換するには適切な個数のグリシン（Ｇ）・オリゴマー鎖で連結する方法が自然と想起されるものとして知られており（特開２００２−１１２７８２号公報参照。）、計算でもそれと同じように扱うことができるが、この方法を常に採用できるとは限らない。即ち、例えば最初の鎖のＣ末端及び／又は次の鎖のＮ末端が複合体界面に位置している場合はグリシン・オリゴマーを割り込ませることができない。また、ＮＭＲの距離情報に基づくディスタンス・ジオメトリー法と同様な方法で複合体を形成させることも可能であるが、この方法は、多くのインプット・データを要することになるので計算準備が煩雑である。そこで、簡便な構造構築方法が求められている。
本発明が解決しようとする課題は、上記の状況に鑑みて、医薬品等の開発に特に重要な鍵となる、複数鎖を有する任意のタンパク質の立体構造を精度良く、かつ従来法に比べてはるかに効率的に構築する方法を提供することにある。また、リガンド分子の多様な改変や受容体等のタンパク質の改変を迅速かつ効率的に行う方法を提供することにもある。更に、本発明は複数鎖から成るタンパク質モデルの構築により、タンパク質−タンパク質間の相互作用様式を解明し、それ等相互作用の認識機構を明確化することで、遺伝病の原因を特定したり、それ等に関連する薬物の開発を促進すること等を目的とする。
２．課題を解決するための手段
本発明者等は、複数鎖を有するタンパク質の任意のアミノ酸配列が与えられたとき適切な立体構造モデルを得るべく検討を重ねた結果、下記（１）−（１０）の方法及びそのためのコンピュータプログラムを見出し、或いは開発した。
（１）任意のアミノ酸配列が与えられた場合において、従来の経験的なモデリング、例えばホモロジーモデリング法の適用範囲を単数鎖（単鎖）のタンパク質から複数鎖より成るタンパク質に拡張（拡張型モデリング法）すべく、コンピュータプログラムにおける当該複数鎖の入力ファイル形式を擬似単鎖（単数鎖）としての形態を有するよう修正し、一方ポテンシャル計算式では複数鎖であることを考慮に入れて、そのタンパク質主鎖及び側鎖の３次元座標から成る立体構造モデルを単数鎖と同様に全自動又は手動的に構築する方法及びそのためのコンピュータプログラム。
目的タンパク質を構成する鎖の数と参照タンパク質を構成する鎖の数は同一であることが望ましい。例えば、目的タンパク質を構成する鎖の数が２個の場合、この立体構造構築のため選択されるされる参照タンパク質（候補）には２個の鎖を有する参照タンパク質を採用することができる。また、同様に目的タンパク質を構成する鎖の数が３個の場合、この立体構造構築のために選択される参照タンパク質（候補）には３個の鎖を有する参照タンパク質を採用することができる。しかしながら必ずしも、両者（参照タンパク質と目的タンパク質）の鎖数が同一である必要はなく、目的タンパク質の立体構造を構築するタンパク質或いはその部分（複数鎖）を参照タンパク質に包含しておれば、それを参照タンパク質として採用することも可能である。
それぞれ複数鎖のタンパク質を１本鎖であるようにコンピュータの配列のファイル形式を修正する（擬似単数鎖）が、擬似単数鎖にする方法には特に制限は無い。例えば、複数鎖より成る目的タンパク質と参照タンパク質とのアミノ酸配列アライメントを、各タンパク質鎖（例えば、ポリペプチド鎖）の境界に識別記号を入れることで、擬似単数鎖としてのアライメントに修正する方法は、当該プログラムを実行する上で簡便かつ有益である。この他にも擬似１本鎖に変換し得る方法が存在する。例えば、当該各タンパク質鎖の境界を、順次、残基番号で登録したファイルを別途作成する等して、コンピュータソフトプログラムへ変数として渡す方法によっても実現可能であるが、両者は各タンパク質鎖の境界を指し示すために用いたファイル形式が異なるだけで、内容自身は全く同等である。同様に、配列アライメント操作自体をコンピュータソフトプログラム内部で行い、その結果から得られる各タンパク質鎖の境界位置又は境界を示す識別数字等を内部変数として利用する方法も好適である。これ等の方法を、本発明において当然用いることができる。従って、これ等の内容も上記開発内容（本発明の内容）に含まれる。
また、当該複数鎖の入力ファイル形式を擬似単数鎖としての形態を有するよう修正することにより、通常、鎖数の増加によって複雑化する入力ファイル形式とタンパク質鎖間の相互作用を、常に最も単純な入力形式と最も明解なポテンシャル計算式で記述することができる。
（２）上記方法において、配列アライメントの修正に際して、好ましくは各種既存のコンピュータソフトウエアであるＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ、ＬＩＢＲＡ、ＦＡＭＳ、ＲＢＳ−ＢＬＡＳＴ、ＩＭＰＡＬＡ、ＣｌｕｓｔａｌＷ、ＨＭＭＥＲ、及びＢＩＯＣＥＳ等のアライメント出力結果を用いて、複数鎖から成る参照タンパク質のアミノ酸配列を選択し、配列アライメントの各アミノ酸配列末尾（Ｃ末端）の次にアミノ酸記号以外の区切り文字（例えば、「Ｕ」等）を付加することで、モデリング時に各鎖間の境界を考慮し、あたかも単数鎖として取り扱うことにより、鎖数の拡張性を有した配列アライメントに修正する方法。尚、前記したように、区切り文字の付加は、複数の鎖を全て相互に結合して単数鎖として扱えるようにするものであり、単数鎖にしたときに結合する最後の鎖のＣ末端部（単数鎖にしたときの最末端部）には必要がない。また、付加しても最末端部の区切り文字を無視することができる。
このようなフォーマットを持つ配列アライメントを準備することにより、タンパク質の鎖数の上限（現在はデフォルトで３６個まで）を、原理的には（計算機上のメモリー制限内で）無限に拡張することができる。
（３）拡張型ホモロジーモデリング等、拡張型モデリング法において、ポテンシャル計算を行う際に、上記（２）の方法で修正した配列アライメントから各タンパク質鎖（例えば、ポリペプチド鎖）のＣ末端残基番号を判定し、その境界における化学結合ポテンシャル及び化学結合角ポテンシャルの切断、並びに当該境界における原子間相互作用ポテンシャルの付加を行うことにより、Ｃα原子座標の構築と最適化を全鎖について順次又は同時に行う、又は行うことができる、モデル構造を構築する方法及びそのためのコンピュータプログラム。
（４）拡張型ホモロジーモデリング等、拡張型モデリング法において、ポテンシャル計算を行う際に、上記（２）の方法で修正した配列アライメントから各タンパク質鎖（例えば、ポリペプチド鎖）のＣ末端残基番号を判定し、その各タンパク質主鎖の境界における化学結合ポテンシャル、化学結合角ポテンシャル、及び化学結合ねじれ角ポテンシャルの切断、並びに当該境界における原子間相互作用ポテンシャルの付加を行うことにより、Ｎ、Ｃα、Ｃ及びＯ（カルボニル酸素）の主鎖原子座標とアミノ酸残基側鎖座標の構築と最適化を全鎖について順次又は同時に行う、又は行うことができる、モデル構造を構築する方法及びそのためのコンピュータプログラム。
尚、図１に当該タンパク質の全アミノ酸配列と各タンパク質鎖（例えば、ポリペプチド鎖）におけるＣ末端残基の通し番号（以下、「ｋ_Ｎ」と称する。Ｎはタンパク質鎖の番号）の関係を、図２〜８には上記（３）及び（４）の方法の具体的な説明を示した。即ち、Ｃα原子座標の構築と最適化の各プロセスにおいては、図２〜４に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間の化学結合ポテンシャル、化学結合角ポテンシャルの切断、及び原子間相互作用ポテンシャルの付加を、Ｎ＝１から全鎖数−１（Ｍ−１）までの全タンパク質鎖について順次又は同時に行い、主鎖原子座標の構築と最適化の各プロセスにおいては、図５〜８に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間の化学結合ポテンシャル、化学結合角ポテンシャル、化学結合ねじれ角ポテンシャルの切断、及び原子間相互作用ポテンシャルの付加を、Ｎ＝１から全鎖数−１（Ｍ−１）までの全タンパク質鎖について順次又は同時に行うことにより、モデル構造の妥当性を向上させることができ、その上計算過程の単純化、効率化を実現することができる。
（５）拡張型モデリング法、特に拡張型ホモロジーモデリング法或いはスレッディング法を用いる際に、タンパク質−タンパク質相互作用部位のアミノ酸残基の原子座標を正確に決定することにより、タンパク質−タンパク質相互作用モデルを構築する方法。
（６）拡張型ホモロジーモデリング等、拡張型モデリング法を用いる際に、複数鎖を構成する或るタンパク質コンポーネント（成分）がペプチド性のリガンド（アミノ酸誘導体、ペプチド誘導体等含む。）を含む場合に、このリガンドのアミノ酸配列を多様に改変したデータセットを作成することにより、目的タンパク質に結合し得る多様なペプチド性のリガンドモデルを構築する方法。
（７）拡張型モデリング法、特に拡張型ホモロジーモデリング法を用いる際に、複数鎖を構成するコンポーネント（成分）がペプチド性のリガンド（アミノ酸誘導体、ペプチド誘導体等含む。）を含む場合に、このリガンドのアミノ酸配列を固定し、周囲のタンパク質鎖のアミノ酸配列を多様に改変したデータセットを作成することで、特定のリガンドに結合し得る多様なリガンド受容体（タンパク質）モデルを構築する方法。
（８）拡張型モデリング法、特に拡張型ホモロジーモデリング法を用いる際に、複数鎖を構成するコンポーネント（成分）が通常のタンパク質（鎖）である場合に、タンパク質同士の認識部位のアミノ酸配列を多様に改変したデータセットを作成して、タンパク質の機能効率を上昇させたり、低下させたりすることにより、有用なタンパク質モデルを構築する方法。
（９）単数鎖の場合において、タンパク質を領域や部品に切り分けることにより複数鎖とみなした後、上記（１）〜（８）記載事項を適用して擬似単数鎖に復元することを可能にすることにより、立体構造モデルの高精度化を図る方法。
（１０）上記（１）〜（９）によって構築されるタンパク質の立体構造モデル、リガンドモデル、及びリガンド受容体（タンパク質）モデルから成るデータベースに対して、下記内容を閲覧、検索できるコンピュータソフトプログラムとそのプログラムを載せたコンピュータ。
閲覧及び検索の対象は、
イ．複数鎖より成る目的タンパク質の遺伝子同定コード又はタンパク質同定コード、１行程度の機能説明、目的アミノ酸配列、及び上記（１）〜（９）により構築される３次元立体構造（の座標）；
ロ．参照タンパク質の遺伝子同定コード又はタンパク質同定コード、１行程度の機能説明、参照アミノ酸配列、及び参照タンパク質の３次元立体構造（の座標）；並びに
ハ．目的配列と参照配列間のアライメント結果、ホモロジー値、及びＥ−ｖａｌｕｅを記述してまとめたデータベース構造。
即ち、本発明は、一つの形態として、
与えられた任意のアミノ酸配列を有する単数鎖より成るタンパク質の立体構造を構築する経験的なモデリング法を拡張して（拡張型モデリング法）、与えられた任意のアミノ酸配列を有する複数鎖より成るタンパク質（目的タンパク質）の立体構造を構築する方法であって、コンピュータソフトプログラムにおける当該複数鎖の入力ファイル形式を擬似単数鎖としての形態を有するように修正し（配列アライメントの修正）、コンピュータソフトプログラムにおけるポテンシャル計算式の計算では複数鎖を有するものとして当該モデリング法によりその立体構造を構築することに特徴を有する複数鎖を有するタンパク質の立体構造構築方法に存する。尚、本発明方法における擬似単数鎖としての入力ファイル形式とポテンシャル計算式での複数鎖としての取り扱いはＦＡＭＳ以外の経験的モデリング法や、更には非経験的モデリング法にも適用可能である。
当該経験的なモデリング法としては、ホモロジーモデリング法（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）及び／又はスレッディング法を挙げることができる。
上記方法においては、全自動的に、或いは手動的に目的とする立体構造を構築することができる。
上記方法において当該修正の内容として、構造が既知で目的配列と同一本数の複数のタンパク質鎖（ポリペプチド鎖、リガンド等）を含む配列を選択し、これを擬似単数鎖として扱うことを挙げることができる。例えば、目的タンパク質がポリペプチド鎖で構成される場合、これと同数のポリペプチド鎖から成る参照タンパク質のアミノ酸配列を選択し、これより、タンパク質を構成する複数のポリペプチド鎖（アミノ酸配列）を一の鎖のＮ末端部（頭部）と他の鎖のＣ末端部（尾部）で順に結合した形での単数鎖（擬似単数鎖）として扱うことを含む。このとき、目的配列を構成する複数のポリペプチド鎖についても、同様に順に結合した形での単数鎖として扱うことができる。例えば、複数鎖を構成する各アミノ酸配列のＣ末端部（末尾）に区切り文字、好ましくはアミノ酸記号以外の区切り文字（例えば、文字「Ｕ」）を付加することにより、例えばコンピュータ上で、単数鎖のタンパク質として取り扱うことができるようにすること（配列アライメントの修正）が好ましい。この場合、順次結合する複数の鎖のうち、最後の鎖、即ちタンパク質がＮ本のタンパク質鎖で構成されている場合Ｎ番目のタンパク質鎖については更に鎖を結合する必要が無いので、このＮ番目のタンパク質鎖（最後の鎖）のＣ末端部には前記区切り文字を付加する必要がない。また、これを鎖の最末端部に付加したとしてもコンピュータ上等でこれを無視することもできる。このように、モデリング時に各鎖間の境界を考慮し、あたかも単数鎖として取り扱うことにより、鎖数の拡張性を有した配列アライメントを利用することができる。
尚、目的タンパク質や参照タンパク質を構成するタンパク質鎖には、実質的にタンパク質を構成する鎖であると認められ、又はタンパク質を構成することができる鎖であれば特に制限は無い。その代表例としてポリペプチド鎖を挙げることができる。「ポリペプチド鎖」とは、アミノ酸及び／又はその誘導体（塩、エステル誘導体等）によるペプチド結合（アミド結合）のみで構成されるものに限らずポリペプチド残基を主として含む鎖であればよい。その他に、タンパク質鎖として、アミノ酸誘導体（例えば、天然若しくは非天然アミノ酸及びその誘導体）、ペプチド誘導体、薬剤成分、核酸類、糖類、有機金属化合物等有機化合物や、金属酸化物及びそのイオン、金属及びそのイオン等無機化合物から選択される物質を挙げることができる。これ等の場合、リガンドとして存在し、又は選択されることが多い。
本発明では、目的タンパク質又は参照タンパク質として使用するタンパク質は複数の鎖で構成され、このような複数のタンパク質鎖には、前記のようなポリペプチド鎖が一つ又は複数含まれ、タンパク質を形成する。例えば、前記複数のポリペプチド鎖のみで構成されるタンパク質、一つ又は複数のポリペプチド鎖と例えば前記リガンドとして存在するアミノ酸誘導体、ペプチド誘導体等の物質を鎖として含むタンパク質を例示することができる。
当該複数鎖の中に、ポリペプチド鎖を複数含むタンパク質の場合、これ等の複数のポリペプチド鎖は、ヘテロ鎖でもホモ鎖でもよい。即ち、当該複数のポリペプチド鎖のアミノ酸配列が全く同一である（ホモ鎖）ような複数のポリペプチド鎖を含んでいてもよいし、相互にそのアミノ酸配列を異にする（ヘテロ鎖）複数のポリペプチド鎖を含んでいてもよい。
上記方法において、参照タンパク質が立体構造データベースから検索され、参照タンパク質のアミノ酸配列と目的配列との配列アライメントを行う方法を含む。
参照タンパク質を検索しアライメントを出力するソフトウエアとして、ＦＡＭＳ、ＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ、ＬＩＢＲＡ、ＲＢＳ−ＢＬＡＳＴ、ＩＭＰＡＬＡ、ＣｌｕｓｔａｌＷ、ＨＭＭＥＲ、ＢＩＯＣＥＳ（Ｐｒｏｔｅｉｎ Ｅｎｇｉｎｅｅｒｉｎｇ，ｖｏｌ．２，Ｎｏ．５，ｐｐ３４７−３５１，１９８９参照。）等を好ましいものとして採用することができる。
上記方法において、修正された配列アライメントとして、各タンパク質鎖（例えば、ポリペプチド鎖）のアミノ酸配列のＣ末端部に区切り文字を有するように配置した書式を有し、区切り文字で区切られたアライメント区間ごとに参照タンパク質ＩＤを指定するファイルを使用することで、同種又は異種の参照タンパク質のアミノ酸配列を記入したマルチプルアライメントに対処可能にして、それ等の重ね合わせにより平均構造を示せるようにするのが好ましい。
上記方法において、修正後の当該擬似単数鎖から各タンパク質鎖（例えば、ポリペプチド鎖）の末端残基番号を判定し、その境界における化学結合角ポテンシャルの切断、及び当該境界における原子間相互作用ポテンシャルの付加を行うことで、立体構造データベース、好ましくはＰＤＢ等、及び／又は当該立体構造データベースを類似構造の重複を回避するように加工し、又は整理したデータベースに基づきＣα原子座標及び主鎖原子座標の構築を行い、シミュレーティッドアニーリング法、分子力学計算、及びモンテカルロ法の少なくとも１種により疑似エネルギー値を表す目的関数の最小化（最適化）を行うことが好ましい。この結果、Ｃα原子座標の構築と最適化を全鎖について順次又は同時に行い、目的とするモデル構造を構築することができる。
例えば、ポテンシャル計算を行う際に、上記の方法で修正した配列アライメントから各タンパク質鎖（例えば、ポリペプチド鎖）のＣ末端残基番号を判定し、その各タンパク質主鎖の境界における化学結合ポテンシャル、化学結合角ポテンシャル、及び化学結合ねじれ角ポテンシャルの切断、並びに当該境界における原子間相互作用ポテンシャルの付加を行うことにより、Ｎ、Ｃα、Ｃ及びＯ（カルボニル酸素）の主鎖原子座標とアミノ酸残基側鎖座標の構築と最適化を全鎖について順次又は同時に行うことができる。更に、上記モデリング法を用いる際に、タンパク質−タンパク質相互作用部位のアミノ酸残基の原子座標を正確に決定し、タンパク質−タンパク質相互作用モデルを構築することができる。
上記方法において、目的タンパク質を構成する複数鎖のうち少なくとも二つの鎖がポリペプチド鎖等のタンパク質鎖である場合、タンパク質−タンパク質相互の認識部位に位置するアミノ酸残基のそれぞれについて２０種類のアミノ酸の可能な組合せに基いて、ポテンシャルエネルギー値を指標として、適合性等の優れた、又は劣った改変体のデータセットを作成して、当該少なくとも二つのタンパク質鎖の各々のタンパク質の機能を上昇させた、又は低下させた立体構造を構築できるようにすることが好ましい。
目的タンパク質を構成する複数鎖のうち少なくとも一つの鎖がアミノ酸誘導体、例えばβＡｓｐ及びγＧｌｕ等の非天然アミノ酸若しくはその誘導体、又はペプチド誘導体（ペプチド性リガンド）であり（複数鎖の少なくとも一つの鎖はポリペプチド鎖である。）、かつ参照タンパク質において該当するリガンド分子と類似した化学構造を持つものである場合に、目的タンパク質の当該誘導体を新たな残基名及び１文字コードで定義し、かつ参照タンパク質の当該リガンドを更に新たな残基名及び１文字コードで定義したアライメントを手動的又は自動的に作成して、当該リガンドの配列を構成するそれぞれの残基について２０種類のアミノ酸及びその誘導体の可能な組み合わせに基づいて、ポテンシャルエネルギー値が低い順にランク付けした場合の上位幾つかを、受容体タンパク質の結合部位への適合性の優れた改変体として格納したアミノ酸誘導体又はペプチド誘導体リガンドモデルデータセットを構築できるようにすることが好ましい。このリガンドのアミノ酸配列を多様に改変したデータセットを作成することにより、目的タンパク質に結合し得る多様なペプチド性のリガンドモデルを構築することができる（後述実施例３参照。）。
同様に、目的タンパク質を構成する複数鎖の鎖成分の中にペプチド性リガンドが存在する場合（複数鎖のうち少なくとも一つはポリペプチド鎖である。）、当該リガンドのアミノ酸配列を固定し、当該リガンドを認識する部位に位置するアミノ酸残基のそれぞれについて２０種類のアミノ酸の可能な組み合わせに基づき、ポテンシャルエネルギー値を指標として、上位複数の受容体タンパク質の結合部位への適合性の優れた改変体のデータセットを作成して、当該リガンドに結合し得る多様なリガンド受容体タンパク質の立体構造を構築できるようにすることが好ましい。このようにして、特定のリガンドに結合し得る多様な受容体モデルを構築することができる。
上記方法において、当該複数鎖としては、単数鎖ポリペプチドを領域や部品に切り分けたもので、擬似単数鎖に復元できるようにしたものでもよい。タンパク質を領域や部品に切り分けることにより複数鎖と見なした後、上記モデリング法を適用して擬似単数鎖に復元することにより立体構造モデルの高精度化を図ることができる。
上記方法においては、当該目的タンパク質（又は目的配列）中に、複数鎖成分として（複数鎖を構成する鎖として）通常のアミノ酸でもなく、また通常のアミノ酸が複数結合したペプチドでもない物質であって、好ましくは立体構造データベース（ＰＤＢ等）に登録されている物質を含むことができる。その例として、非天然アミノ酸、薬剤成分、核酸類、糖類、有機金属化合物等有機化合物や、金属酸化物、そのイオン、金属及びそのイオン等無機化合物等のヘテロ成分を挙げることができる。当該タンパク質には、それを構成する複数鎖の少なくとも一つとしてこのような物質を（例えば、リガンド等で）含む場合でもポリペプチド鎖が一つ又は複数含まれて、当該タンパク質を形成している。
上記方法として、立体構造データベースの中から、目的配列に適した参照タンパク質を検索し、検索した複数の参照タンパク質のアミノ酸配列との配列アライメントを行うステップ；目的配列に対してＥ−ｖａｌｕｅの小さい参照タンパク質のアミノ酸配列を選択するステップ；及び参照タンパク質に含まれる鎖のアミノ酸配列のＣ末端部（末尾）に区切り文字を付加し、目的配列の対応する位置にも区切り文字を付加する（配列アライメントの修正）ステップを含むことができる。
尚、前記したように区切り文字の付加は、複数の鎖を全て相互に結合して単数鎖として扱えるようにするものであり、単数鎖にしたときに結合する最後の鎖のＣ末端部（単数鎖にしたときの最末端部）には、更に結合する鎖が存在しないので必要がない。付加の必要はないが同様に付加して、付加された最末端部の区切り文字を無視することもできる。
更に、アライメントの情報に基いて、前記の如く、目的配列のアミノ酸中の構成原子の一つであるＣα原子について、前記参照タンパク質のアミノ酸配列を選択するステップで定められた参照構造から座標を取得し、上記発明方法（例：請求の範囲１１参照。）において使用するＣαの原子座標を最適化するステップ；得られたＣα座標に主鎖原子座標を立体構造データベースより付加し、上記発明方法において使用する主鎖原子座標を最適化するステップ；及び得られた主鎖原子座標に側鎖原子座標を立体構造データベースより付加し、上記発明方法において使用する側鎖原子座標を最適化するステップを含むことができる。
上記ポテンシャル計算式としては、後記表１に示されているポテンシャル計算式を採用することができる。好ましくは、下記の通りである。
全鎖数＝Ｍのときのポテンシャル計算式については、Ｎがポリペプチド鎖等のタンパク質鎖の番号を、ｋ_ＮがＮ番目のタンパク質鎖（ポリペプチド鎖等）におけるＣ末端残基の通し番号をそれぞれ表し、ｉ＝１，．．．，Ｍ−１を簡略化して、ｉ＝１，Ｍ−１と表記するとき、
（Ａ）Ｃα原子座標の構築と最適化プロセスにおける計算では、擬似化学結合ポテンシャルのｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の場合を含めず、擬似化学結合角ポテンシャルのｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}＋１の場合を含めない。また、原子間相互作用ポテンシャルの場合は、ｉ＝ｋ_Ｎ−１ならｊ＞ｉ＋１、ｉ＝ｋ_Ｎならｊ＞ｉを各々追加する。
（Ｂ）主鎖原子座標の構築と最適化プロセスにおける計算では、化学結合ポテンシャルにはｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}のときのＣｉとＮｉ＋１間の結合を含めず、化学結合角ポテンシャルには、Ｃ、Ｏをそれぞれカルボニルの炭素原子と酸素原子、Ｃαをα炭素原子、Ｎを窒素原子とするとき、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の角度Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１、Ｏｉ−Ｃｉ−Ｎ_ｉ＋１及びＣｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１を含めず、更に、化学結合ねじれ角ポテンシャルには、ｉ＝ｋ_{Ｎ（ Ｎ＝１，Ｍ−１）}のときの角度Ｎｉ−Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１、Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１、及びＣｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１−Ｃ_ｉ＋１を含めない。また、原子間相互作用ポテンシャルには、原子間距離をｒとするとき、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の、ｒ_ｉｊ∈｛ｒ_{Ｎｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｎｉ＋１}；ｒ_{Ｃｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｃβｉ＋１}；ｒ_{Ｃｉ，Ｃｉ＋１}；ｒ_{Ｏｉ，Ｎｉ＋１}；ｒ_{Ｏｉ，Ｃαｉ＋１}｝についてｒ_ｉｊ≦規定値の場合を追加する。
即ち、拡張型モデリング法において、ポテンシャル計算を行う際に、前記方法で修正された配列アライメントから各タンパク質鎖（ポリペプチド鎖等）のＣ末端残基番号を判定し、その各タンパク質鎖の主鎖の境界における化学結合ポテンシャル、化学結合角ポテンシャル、及び化学結合ねじれ角ポテンシャルの切断、並びに当該境界における原子間相互作用ポテンシャルの付加を行うことで、Ｎ、Ｃα、Ｃ及びＯ（カルボニル酸素原子）の主鎖原子座標とアミノ酸残基側鎖座標の構築と最適化を全鎖について順次又は同時に行い、モデル構造を構築することができる。
尚、図１には当該タンパク質の全アミノ酸配列と各タンパク質鎖（ポリペプチド鎖等）におけるＣ末端残基の通し番号の関係を示したものである。全鎖数＝Ｍ、各ｋ_Ｎは区切り文字Ｕより識別する。Ｎはタンパク質鎖（ポリペプチド鎖等）の番号を表す。
図２〜８には上記配列アライメントの修正方法の具体的な説明を示した。即ち、Ｃα原子座標の構築と最適化の各プロセスにおいては、図２〜４に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間の化学結合ポテンシャル、化学結合角ポテンシャルの切断、及び原子間相互作用ポテンシャルの付加を、Ｎ＝１から全鎖数−１（Ｍ−１）までの全タンパク質鎖（ポリペプチド鎖等）について順次又は同時に行い、主鎖原子座標の構築と最適化の各プロセスにおいては、図５〜８に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間の化学結合ポテンシャル、化学結合角ポテンシャル、化学結合ねじれ角ポテンシャルの切断、及び原子間相互作用ポテンシャルの付加を、Ｎ＝１から全鎖数−１（Ｍ−１）までの全タンパク質鎖について順次又は同時に行うことにより、モデル構造の妥当性を向上させることができることに加え、計算過程の単純化及び効率化を実現することができる。
入力ファイル形式としては特に制限は無いが、例えば、目的タンパク質を参照タンパク質からモデリングする場合の入力ファイル形式の一般式として、下記内容を採用することができる。
１行目は＞の後に、目的タンパク質のＩＤを記載する。２行目は目的タンパク質のアミノ酸配列を各タンパク質鎖（ポリペプチド鎖等）の末端残基の次に区切り文字（列）を加えて、全てのタンパク質鎖について改行無しに記載する。３行目は＞の後に、参照タンパク質のＩＤを記載する。４行目は参照タンパク質のアミノ酸配列を各タンパク質鎖（ポリペプチド鎖等）の末端残基の次に区切り文字（列）を加えて、全てのタンパク質鎖について改行無しに記載する。２行目と４行目のアミノ酸配列の並べ方については、前記の如く得られるアライメントを利用し、同一番目の区切り文字（列）の位置は２行目と４行目で常に一致させるのが好ましい。
これを例示すると、以下の通りである。
＞（目的タンパク質のＩＤ）；
（第一鎖の配列）（区切り文字）（第二鎖の配列）（区切り文字）…（区切り文字）（第Ｍ鎖の配列）；
＞（参照タンパク質のＩＤ）；
（第一鎖の配列）（区切り文字）（第二鎖の配列）（区切り文字）…（区切り文字）（第Ｍ鎖の配列）。
尚、上記括弧は見易くするためのものであり、実際には記載せず、連続した文字列とするのが好ましい。
上記本発明の方法において、前述したように、複数鎖より成る目的タンパク質と参照タンパク質とのアミノ酸配列アライメントを、各タンパク質鎖（ポリペプチド鎖等）の境界に識別記号を入れることで、擬似単数鎖としてのアライメントに修正する方法は、当該プログラムを実行する上で簡便かつ有益であるが、他に変え得る方法が無い訳ではなく、例えば、各タンパク質鎖の境界を、順次、残基番号で登録したファイルを別途作成する等して、コンピュータソフトプログラムへ変数として渡す方法によっても実現可能であるが、両者は各タンパク質鎖の境界を指し示すために用いたファイル形式が異なるだけで、内容自身は全く同等のものである。同様に、配列アライメント操作自体をコンピュータソフトプログラム内部で行い、その結果から得られる各タンパク質鎖（ポリペプチド鎖等）の境界位置又は境界を示す識別数字等を内部変数として利用する方法も全く同等のものであり、これ等の方法を、本発明として当然用いることができる。
本発明は、与えられた任意のアミノ酸配列を有する単数鎖より成るタンパク質の立体構造を構築する経験的なモデリング法を拡張して（拡張型モデリング法）、与えられた任意のアミノ酸配列を有する複数鎖より成るタンパク質（目的タンパク質）の立体構造を構築する方法であって、当該目的タンパク質及び選択された参照タンパク質について、それぞれその中に含まれ、相互に対応する複数鎖の各アミノ酸配列を、Ｎ末端部とＣ末端部とで順次結合した状態の単数鎖とみなした後、このように得られた擬似単数鎖の参照配列と擬似単数鎖の目的配列との配列アライメントを行い対応関係の確認を行った後、目的配列中のアミノ酸残基中の構成原子の一つであるＣα原子を配置し、その間をアミド結合でつなぎ、更に側鎖を付加してその他の構成原子についても座標を構築し、最適化を行い、当該モデリング法によりその立体構造を構築することに特徴を有する複数鎖を有するタンパク質の立体構造構築方法にも存する。
選択される参照タンパク質としては、目的タンパク質に対して、ポリペプチド鎖等のタンパク質鎖のアミノ酸配列が類似したもの、特に統計的に有意に類似したものを好ましく選択することができる。
また、前記本発明の立体構造構築方法において、目的タンパク質の立体構造構築を、得られたアライメントの情報に基いて、目的タンパク質主鎖アミノ酸中のＣα原子について選択された参照タンパク質の立体構造から座標を取得し、目的関数を最小化するようにＣαの原子座標を最適化し、最適化されたＣαの原子座標に主鎖の他の原子座標（Ｃβの原子座標を含む。）を付加して目的関数を最小化するように主鎖の原子座標を最適化し、最適化された主鎖の原子座標に側鎖の他の原子座標を付加し目的関数を最小化するように側鎖の原子座標を最適化することにより行うことができる。
本発明は、別の形態として、上記本発明方法で構築されたことに特徴を有するタンパク質の立体構造モデルにも存する。更には、前記発明方法の何れかの方法を含むことに特徴を有するプログラム（ｎｅｗＦＡＭＳ）、又は当該プログラムを載せたことに特徴を有するコンピュータにも存する。
本発明は、更に別の形態として、上記本発明方法で構築されたタンパク質の立体構造モデル、並びに上記本発明方法の中において使用されるリガンドモデル及び同様に使用されるリガンド受容体タンパク質の立体構造モデルから成るデータを固定化し、組み合わせたことに特徴を有する拡張型モデリング法に使用可能なデータベースにも存する。
本発明は、更なる別の形態として、上記本発明方法で構築されたタンパク質の立体構造モデル、並びに上記本発明方法の中において使用されるリガンドモデル及び同様に使用されるリガンド受容体タンパク質の立体構造モデルのデータを、コンピュータで閲覧又は検索できるようにデータベース化したことに特徴を有するデータベースにも存する。
本発明は、もう一つ別の形態として、下記内容をコンピュータで閲覧又は検索できるようにしたことに特徴を有するデータベース構造にも存する：
複数鎖より成る目的タンパク質の遺伝子同定コード又はタンパク質同定コード、１行程度の機能説明、目的アミノ酸配列、及び当該目的タンパク質の３次元立体構造（の座標）；
参照タンパク質の遺伝子同定コード又はタンパク質同定コード、１行程度の機能説明、参照アミノ酸配列、及び当該参照タンパク質の３次元立体構造（の座標）；並びに
目的配列と参照配列間のアライメント結果、ホモロジー値、及びＥ−ｖａｌｕｅ。
本発明は、更に別の形態として、上記何れかの、データベースの内容を閲覧又は検索でき、又はデータベース構造を利用できるようにしたことに特徴を有するコンピュータソフトプログラム又は当該プログラムを搭載したコンピュータにも存する。
本発明は、更に、前記本発明の何れかの方法で構築された立体構造データベースの中から、閲覧したいタンパク質を、生物種固有の任意の記号、タンパク質コード名、参照タンパク質名、１行程度の機能説明の文字列との部分一致等の、論理積検索によって、目的とするタンパク質にたどり着けるように、好ましくは簡便に事前の知識無く目的とするタンパク質にたどり着けるように設計されたことに特徴を有するインターフェースにも存する。
上記方法により、複数鎖を有するタンパク質の立体構造を規定する原子座標が提供される。ここで、本明細書においてＣα原子とは、各アミノ酸残基の骨格の中心となる炭素原子を意味し、グリシンを除く各アミノ酸残基のＣα原子は光学活性の性質をもつ。Ｃβ原子とは、Ｃα原子の側鎖側に結合する炭素原子を意味する。また、Ｃ原子とは、Ｃα原子に結合するカルボニル基の炭素原子を意味する。
表１で、Ｎはポリペプチド鎖等タンパク質鎖の番号、ｋ_ＮはＮ番目のタンパク質鎖（ポリペプチド鎖等）におけるＣ末端残基の通し番号を表し、全鎖数＝Ｍのとき、ｉ＝１，．．．，Ｍ−１を、ｉ＝１，Ｍ−１と表記する。図１には当該タンパク質の全アミノ酸配列とｋ_Ｎの関係を示した。また、表１中の各ポテンシャル計算式について、各タンパク質鎖（ポリペプチド鎖等）の境界における扱いについては、図２〜８に具体的に説明を示し、また次節において詳細に説明する。
［表１］
各タンパク質鎖の境界におけるポテンシャル関数の計算条件
従来法（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００等参照。）による。この中で、従来型の単数鎖モデリングでは知られておらず、また示唆されてもいない内容で、本発明の複数鎖タンパク質のための拡張型モデリング法において新規に改訂された部分が、下記ポテンシャル式の「但し、．．．」（但し書き）によって示されている。即ち、下記但し書き以外の部分は従来法のモデリング法が用いられる。従って、上記但し書きの内容については、本発明で使用する擬似単数鎖モデリングにおいて新規に採用することができる修正内容に該当する。尚、特別の説明が無い限り、計算式中の文字の意味は、全て上記従来法の文献に従い、計算式中の定数は、合理的且つ妥当な範囲で変動し得るものであり、下記具体的値に限定されるものではない。ここで、当該従来法の文献の全内容は参考により本明細書中にその一部として組み込まれる。
（Ａ）Ｃα原子座標の構築と最適化プロセス

（Ｂ）主鎖原子座標の構築と最適化プロセス

前記したように、本発明方法において、当該複数鎖の入力ファイル形式とポテンシャル計算式をそれぞれ擬似単数鎖としての形態を有するよう修正することにより、通常、鎖数の増加によって複雑化する入力ファイル形式とタンパク質鎖（ポリペプチド鎖等）間の相互作用を、常に最も単純な入力形式と最も明解なポテンシャル計算式で記述することができる。
発明の実施の形態
以下、本発明の実施の形態について説明する。好適かつ代表的な例としてホモロジーモデリング法を中心に、本発明を詳細に説明するが、代表的な例について説明するものであり、本発明はこれに限定されるものではない。
本発明は、北里大学薬学部生物分子設計学教室（梅山秀明教授）で開発した蛋白質自動モデリングシステム（ＦＡＭＳ：Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。この全内容は参考により本明細書中にその一部として組み込まれる。）の適用範囲を拡張し、改良するものでもある。
本明細書において幾つかの用語を使用するが、特に明記しない限り、次の意味を有する。
「目的タンパク質」とは、Ｘ線結晶解析やＮＭＲ解析等により完全な立体構造が決定されておらず、本発明において立体構造構築の対象とするタンパク質を意味する。このタンパク質のアミノ酸配列を「目的配列」又は「目的アミノ酸配列」と称することがある。この目的タンパク質には、部分構造は解析されているが完全な立体構造が得られていないもの、また、既に機能が特定されているもの、機能が推定されているもの、アミノ酸配列は決定されているが機能は全く不明のもの等も含まれる。「参照タンパク質」とは、その立体構造の詳細がＸ線結晶解析やＮＭＲ解析により既に決定されており、アライメントや原子座標の最適化のために参照するタンパク質を意味する。このタンパク質のアミノ酸配列を「参照配列」又は「参照アミノ酸配列」と称することがある。
「アライメント」とは、２種類以上のタンパク質のアミノ酸配列があった場合、相互の対応関係をつけることを意味し、その方法は以下の各ステップの説明において詳述する。
「原子座標」とは、３次元座標上で立体構造を記述するものである。それは空間上のある点を原点とする互いに垂直な三方向の相対的な距離であり、タンパク質中に存在する水素原子を除く原子一つ当たりに３個の数字から成るベクトル量である。
図９は、本発明による複数鎖から成るタンパク質の立体構造構築法の一例を示すフローチャートである。
図９に示す通り、この方法（一例）においては、先ず、ステップ１０において、立体構造未知のタンパク質のアミノ酸配列（以下、「目的配列」と称する。）を用意する。ステップ２０において、ＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ又はＬＩＢＲＡを用いて、立体構造データベースの中から参照タンパク質（参照アミノ酸配列）を選定する。また、アミノ酸残基を表す２０種の文字を指標として、目的配列と選定された参照アミノ酸配列とのアライメント（並置）を行う。ステップ３０において、検索結果から一つ又は複数の参照タンパク質のアミノ酸配列を選択する。ステップ４０において、複数鎖の各アミノ酸配列のＣ末端部（末尾）に「Ｕ」という文字を区切りに加入し、同じく目的配列の対応する位置の各アミノ酸配列末尾に順次、「Ｕ」を区切り文字として加入することで配列アライメントを修正する。尚、このとき前記の如く最終鎖に該当する鎖の末端部（擬似単数鎖を構成する鎖の末端部；擬似単数鎖の最末端部）には上記区切り文字を付加する必要は無い。また、これを付加しても最末端部に付加された区切り文字を無視することができる。ステップ５０において、アライメントの情報に基づいて、アミノ酸残基中の構成原子の一つであるＣα原子についてステップ３０において定められた参照構造から座標を取得し、シミュレーティッドアニーリング法によって、後述の種々のポテンシャル項の和からなる目的関数（Ｅ_Ｃα）を最小化するようにＣαの原子座標を最適化する。ステップ６０において、ステップ５０で得られたＣα座標に主鎖原子座標をデータベースより付加し、シミュレーティッドアニーリング法によって、後述の種々のポテンシャル項の和からなる目的関数（Ｅ_ｍａｉｎ）を最小化するように主鎖原子座標を最適化する。ステップ７０において、ステップ６０で得られた主鎖原子座標に側鎖座標をデータベースより付加し、ステップ６０と同様のシミュレーティッドアニーリング法によって立体構造の構築を行う。ステップ８０において、でき上がったモデルの立体構造上の妥当性を後述のように検査し、ステップ９０で最終構造に至る。モデル構造の妥当性の検査は、ループの挿入や欠損部位を除いた領域で、モデルと参照タンパク質の主鎖原子座標の重ね合わせを行い、ｒｍｓｄが、例えば１Å以下であれば、妥当なモデルとした。
本発明では、従来は単数鎖の参照構造を用いた構造最適化しかできなかったステップ５０、６０及び７０を修正し、複数鎖にも対応できるようにしたことにより、従来のＦＡＭＳと異なって発展している。以下、好ましい形態としての各ステップについて更に詳細に説明する。
ステップ１０：構造未知のタンパク質のアミノ酸配列
先ず、構造未知の目的タンパク質のアミノ酸配列（目的配列）を用意する。用いられる目的タンパク質のアミノ酸配列としては、データベースに登録されているもの、配列が初めて解析されたもの等如何なる由来の配列であってもよい。また、部分構造のみが解析されているタンパク質の配列も完全な立体構造の情報を得るために、本発明における立体構造の構築対象とすることができる。用いられるデータベースとしては、例えば、ＧｅｎｅＢａｎｋ：ｆｔｐ：／／ｎｃｂｉ．ｎｌｍ．ｎｉｈ．ｇｏｖ／ｇｅｎｂａｎｋ／ｇｅｎｏｍｅｓ／、ＰＩＲ：ｈｔｔｐ：／／ｗｗｗ−ｎｂｒｆ．ｇｅｏｒｇｅｔｏｗｎ．ｅｄｕ／ｐｉｒ／（Ｎａｔｉｏｎａｌ Ｂｉｏｍｅｄｉｃａｌ Ｒｅｓｅａｒｃｈ Ｆｏｕｎｄａｔｉｏｎ（ＮＢＲＦ））、Ｓｗｉｓｓ Ｐｌｏｔ：ｈｔｔｐ：／／ｗｗｗ．ｅｘｐａｓｙ．ｃｈ／ｓｐｒｏｔ／ｓｐｒｏｔ−ｔｏｐ．ｈｔｍｌ（Ｓｗｉｓｓ Ｉｎｓｔｉｔｕｔｅ ｏｆ Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ（ＳＩＢ）、Ｅｕｒｏｐｅａｎ Ｂｉｏｉｎｆｏｍａｔｉｃｓ Ｉｎｓｔｉｔｕｔｅ（ＥＢＩ））、ＴｒＥＭＢＬ（ＵＲＬ及び管理者共にＳｗｉｓｓ Ｐｌｏｔと同じ）、ＴｒＥＭＢＬＮＥＷ（ＵＲＬ及び管理者共にＳｗｉｓｓ Ｐｌｏｔと同じ）、ＤＡＤ：ｆｔｐ：／／ｆｔｐ．ｄｄｂｊ．ｎｉｇ．ａｃ．ｊｐ（日本ＤＮＡデータバンク）等のデータベースに登録されているヒト（Ｈ．ｓａｐｉｅｎｓ）、ショウジョウバエ（Ｄ．ｍｅｌａｎｏｇａｓｔｅｒ）、線虫（Ｃ．ｅｌａｇａｎｓ）、酵母（Ｓ．ｃｅｒｅｖｉｓｉａｅ）、シロイヌナズナ（Ａ．ｔｈａｌｉａｎａ）等を挙げることができる。
これらのデータベースは単なる例示であり、タンパク質のアミノ酸配列が登録されているものであれば如何なるデータベースを用いることもできる。
ステップ２０：ＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ、ＬＩＢＲＡ等のアライメントソフトによるデータベース検索及び配列アライメント
ステップ１０において用意されたアミノ酸配列に対して、アライメントを行うソフトウエアとしては、例えばＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ（Ｐｏｓｉｔｉｏｎ−Ｓｐｅｃｉｆｉｃ Ｉｔｅｒａｔｅｄ ＢＬＡＳＴ）及びＬＩＢＲＡ等を使うのが好ましい。
ＦＡＳＴＡは目的配列と一致度の高い配列を立体構造データベースから探索し、最終的な目的配列と参照タンパク質との一致度をＥ−ｖａｌｕｅとして算出するプログラムである。ＦＡＳＴＡの詳細は″Ｅｆｆｅｃｔｉｖｅ ｐｒｏｔｅｉｎ ｓｅｑｕｅｎｃｅ ｃｏｍｐａｒｉｓｏｎ″Ｐｅａｒｓｏｎ ＷＲ，（１９９６）Ｍｅｔｈｏｄｓ Ｅｎｚｙｍｏｌ；２６６：２２７−５８に記載されている。
ＰＳＩ−ＢＬＡＳＴはプロファイルアライメントを行うようにプログラムされている。ＰＳＩ−ＢＬＡＳＴの詳細は、″Ｍａｔｃｈｉｎｇ ａ ｐｒｏｔｅｉｎ ｓｅｑｕｅｎｃｅ ａｇａｉｎｓｔ ａ ｃｏｌｌｅｃｔｉｏｎ ｏｆ ＰＳＩ−ＢＬＡＳＴ−ｃｏｎｓｔｒｕｃｔｅｄ ｐｏｓｉｔｉｏｎ−ｓｐｅｃｉｆｉｃ ｓｃｏｒｅ ｍａｔｒｉｃｅｓ″Ｓｃｈａｆｆｅｒ ＡＡ，Ｗｏｌｆ ＹＩ，Ｐｏｎｔｉｎｇ ＣＰ，Ｋｏｏｎｉｎ ＥＶ，Ａｒａｖｉｎｄ Ｌ ａｎｄ Ａｌｔｓｃｈｕｌ ＳＦ，Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ １９９９，１２，１０００−１１に記載されている。
プロファイルアライメントを実行するＰＳＩ−ＢＬＡＳＴは配列の類似度を検出することについて、現時点で最高の性能を備えたツールである。ＰＳＩ−ＢＬＡＳＴはＢＬＡＳＴと呼ばれる類似タンパク質を検索し、アライメントを出力する一連のプログラムの１つであるが、最近同等の性能を発揮するものとして、ＲＢＳ−ＢＬＡＳＴとＩＭＰＡＬＡ（Ａ．Ａ．Ｓｃｈａｆｆｅ ｅｔ ａｌ．，Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ，１５（１２），１０００−１０１１，１９９９参照。）がある。このプログラムは、データベース中の統計的に有意なアライメント関係のみからプロファイル情報を引き出し、アミノ酸配列の部位特異的スコア行列（ｐｏｓｉｔｉｏｎ ｓｐｅｃｉｆｉｃ ｓｃｏｒｅ ｍａｔｒｉｘ；アミノ酸配列上で各残基ごとに、あるアミノ酸残基への置換確率を統計的に示した行列。）を作成する。次に、プログラム内部では目的タンパク質の配列の替わりに作成された部位特異的スコア行列と一致度の高い配列をデータベースから探索し、その部位特異的スコア行列はＥ−ｖａｌｕｅの小ささを限界値として見たときに、有意なアライメントが検出されなくなるまで回数を重ねるごとに逐次更新されてゆく。そして、最終的な部位特異的スコア行列と参照タンパク質との一致度をＥ−ｖａｌｕｅとして算出する。
「Ｅ−ｖａｌｕｅ」とは配列間の適合の間に存在するランダムなバックグラウンドノイズを定量的に記述したものである。それは二つの配列がどれだけマッチしているかを示すものでもあり、配列類似度を表現するスコアに対し指数関数的に減少する性質を持ち、結果を評価する際の閾値を設定する方法として有用である。ＰＳＩ−ＢＬＡＳＴにおいては、Ｅ−ｖａｌｕｅの大きさについて、通常０．１以下、好ましくは０．００１以下の値を持つときがそれに相当する（Ａ．Ａ．Ｓｃｈａｆｆｅ ｅｔ ａｌ．，Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ，１５（１２），１０００−１０１１，１９９９参照。）。
これ等により立体構造データベースＰＤＢから参照タンパク質配列をホモロジー検索し、検索された参照配列と目的配列との配列アライメントを行う。
ここで、「参照タンパク質」とは立体構造データベースから取得される配列と３次元原子座標のデータであり、プロテインデータバンク（ＰＤＢ）として登録された公共データベースから得ることができるものである。
ＰＤＢのデータベースには、２００１年１１月現在２６２４３個の立体構造が登録されているが、一例として９５％以上の配列の一致度があるものは同じカテゴリーと判断し、カテゴリー中からもっとも長い配列のものを、更に同じ長さの場合はＸ線の分解能の最も高い構造をそのカテゴリーの代表とした。本発明で用いた立体構造データベースはこの代表を集めたデータベースである。本発明がなされた現在においては、３９２２個の代表構造をＰＤＢデータベースとして用いている。
ステップ３０：検索結果から、一つ又は複数の参照タンパク質のアミノ酸配列を選択
ホモロジー検索の結果から、目的配列に対して統計的に有意に類似した一つ又は複数の参照タンパク質のアミノ酸配列を選択する。
ステップ４０：複数鎖の各アミノ酸配列末尾に区切り文字Ｕを挿入し、配列アライメントを修正
ステップ３０により、複数鎖より成るタンパク質の目的配列が、複数のタンパク質の参照配列とそれぞれアライメントされた場合、従来は目的配列の各鎖ごとにモデリングを行っていた。しかし、本発明では複数鎖の各アミノ酸配列のＣ末端部（末尾）に区切り文字（例えば、文字「Ｕ」）を挿入し、全ての目的タンパク質鎖（ポリペプチド鎖等）を擬似単数鎖として同時にモデリングするためのアライメントを準備する。計算できるタンパク質鎖の上限は、実際に好ましくはＰＤＢに登録された複数鎖タンパク質を参考に３６個としたが、原理的にはモデリングに必要な個数まで、又は用いる計算機能力の限界まで拡張可能である。
ステップ５０：Ｃα原子の初期座標の構築
ステップ４０からの区切り文字を含んだアライメントを利用して、目的配列と参照配列との比較から、挿入及び欠損のあるアミノ酸残基についての情報を得る。アライメントにおいて連続して三残基以上のアミノ酸残基が対応している領域を選び出し、その領域においては、これらのアミノ酸残基ペアにおいて、目的タンパク質のＣα原子として参照タンパク質と同一のものを当てはめておく。Ｃα原子が求められなかった場合には、ＰＤＢから予め作成した、Ｃα原子から成るペプチド断片のデータベース（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）から座標を当てはめる（図１０参照。）。
ステップ５０（１）：Ｃα原子のシミュレーティッドアニーリング法による最適化
上記ステップ５０で作成されたＣα原子は、前記背景技術の項において説明したシミュレーティッドアニーリングのプロセスを用いて参照タンパク質の座標を参照して得られた目的関数（Ｅ_Ｃα）を用いて最適化される。この目的関数は下記式（１）の通りである。従来法と本発明の重要な違いの一つは、式（１）中の化学結合ポテンシャルＥ_ｌｅｎと化学結合角ポテンシャルＥ_ａｎｇ及び原子間相互作用ポテンシャルであるＥ_ｖｄｗを、ステップ４０で作成されたアライメント中の区切り文字Ｕを参考にして定めたｋ_Ｎの値（Ｎ番目のタンパク質鎖のＣ末端残基の通し番号に等しい。）によって下記のように修正した点である。

Ｅ_ｌｅｎは、配列上隣り合う残基のＣα原子間の距離に関するもので下記式（２）のように設定される。

ここで、Ｄ_{ｉ，ｉ＋１}は残基１と残基ｉ＋１のＣα間距離である。Ｋ_１は定数で例えば、２と設定される。但し、１番目のタンパク質鎖のＣ末端残基ｋ_１と２番目のタンパク質鎖のＮ末端残基ｋ_１＋１の間には化学結合が存在しないので、ｉ＝ｋ_１の場合はＥ_ｌｅｎの計算に含めない。同様に、図２に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_１＋１の間には化学結合が存在しないので、ｉ＝ｋ_Ｎの場合はＥ_ｌｅｎの計算に含めない。（以下、このような操作を「相互作用の切断」と称する。）タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の１番目の区切り文字Ｕ_１からＭ−１番目の区切り文字Ｕ_Ｍ−１までのＵの位置によって指定することができる。
次に、Ｅ_ａｎｇはＣα原子の化学結合角の関数であり下記式（３）の通りである。

但し、ｉ≠ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}，ｉ≠ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}＋１
ここで、θ_ｉ（ｒａｄ）はｉ．ｉ＋１、ｉ＋２番目の残基Ｃα原子の角度である。θ_０はＰＤＢのＸ線構造から（１００／１８０）π（ｒａｄ）と設定される。Ｋ_ａは定数で例えば１とする。但し、結合角のポテンシャルＥ_ａｎｇについても、図３に示したように、Ｅ_ｌｅｎと同様に相互作用の切断操作を行う。即ち、ｉ＝ｋ_Ｎ及びｉ＝ｋ_Ｎ＋１の場合は、Ｅ_ａｎｇの計算に含めない。この相互作用の切断操作は、ｋ_１からｋ_Ｍ−１までのアミノ酸残基について行う。
次にＥ_ｖｄｗはＣα原子間のファンデルワールスポテンシャルで、通常は３残基以上離れているものについて考慮し、下記式（４）の通りである。

但し、ｉ＝ｋ_Ｎ−１ならｊ＞ｉ＋１、ｉ＝ｋ_Ｎならｊ＞１
ここで、Ｄ_ｉ，ｊは、対象となるｉ番目のＣα原子から６Å以内にある原子対ｉ，ｊ間の距離であり、Ｋ_ｖｄｗの値は０．０１（Ｄ_ｉ，ｊ≦３．２Å）又は０．００１（Ｄ_ｉ，ｊ＞３．２Å）と設定した。但し、図４に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間には化学結合が存在しないので、ｉ＝ｋ_Ｎ−１ならｊ＞ｉ＋１、ｉ＝ｋ_Ｎならｊ＞ｉという条件でＥ_ｖｄｗの計算を行わなければならない。タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の区切り文字Ｕ_１からＵ_Ｍ−１までのＵの位置によって指定することができる。
次に、Ｅ_ＳＳは対となってＳ−Ｓ結合を構成するＣｙｓ残基のＣα原子間の距離に関するもので下記式（５）のように設定される。

ここで、Ｄ_ｉ ^ＳＳは、タンパク質鎖内及びタンパク質鎖間でジスルフィド結合を形成するＣｙｓ残基のペア同士のＣα間距離である。本発明では、残基番号ｉには全タンパク質鎖にわたる通し番号を用いるため、タンパク質鎖間のＥ_ＳＳを従来法のポテンシャル関数で扱うことができる。Ｋ_ＳＳは定数で例えば５と設定される。
次に、Ｅ_ｐｏｓはＣα原子の位置に関する関数であり、下記式（６）の通りである。このエネルギー項は、タンパク質のＳＣＲｓ（Ｓｔｒｕｃｔｕｒａｌ Ｃｏｎｓｅｒｖｅｄ Ｒｅｇｉｏｎｓ：後述。）におけるＣα原子の位置を比較的安定に保持させる目的で導入したものである。

ここで、Ｘ_ｉはｉ番目のＣα原子の座標を表し、Ｍ_ｉは構造を基にしたアライメント上で構造的に等価な、即ち３次元座標系で最も近接した位置にあるＣα原子間の平均距離である。残基ｉについてＭ_ｉの値が求められないとき、即ち目的配列のあるアミノ酸残基が参照配列のＣαに対応づけられないときは、Ｍ_ｉの値は１０と設定される。‖・‖の意味する所はノルム（座標ベクトル間の距離）であり、ここで〈ｗ_ｉｘ_ｉ〉は、Ｃα原子の平均座標であり下記式（７）の通りである。

但し、ｊ≠ｉ
ここで、Ｘ^ｊ _ｉはｊ番目の参照タンパク質のｉ番目の残基に対応するＣαの原子座標である。ｗ^ｊ _ｉは、ｊ番目の参照タンパク質のｉ番目のＣα原子の位置に対する重みを表し、Ｗはｊに対するｗ^ｊ _ｉの和である。このｗ^ｊ _ｉは、目的タンパク質の大体の形を決定するため重要なパラメータであるが、これは図１１に示すように、ローカルスペースホモロジー（ＬＳＨ）と呼ばれる着目部位の１２Å以内の空間的近傍の局所的な値によって決定されている（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。ＬＳＨと構造がよく保存されている部位（ＳＣＲｓ：Ｓｔｒｕｃｔｕｒａｌ Ｃｏｎｓｅｒｖｅｄ Ｒｅｇｉｏｎｓ）に存在する残基のペアの比率との相関は図１２に示されているように非常に高い。これは、高いＬＳＨ値を持つときは統計的にＣα原子の位置が参照タンパク質構造と比べて１．０Å以内にあることを意味する。
Ｃα原子は式（１）に従って、シミュレーティッドアニーリング法を用いて、繰り返し最適化される。この最適化の段階でＣα原子の摂動は、例えば１．０Å以内になるように設定する。また、このアニーリングの段階は全てのＣα原子について例えば１００回づつ計算される。そして、温度に相当するパラメータは、２５からスタートし、１回毎に０．５倍にしてゆき、０．０１になるまで減らすようにし、そのパラメータは以後一定とする。
立体構造の重ね合わせの結果、目的配列に対して最も挿入、欠損が少ない参照配列からの構造情報の取得とＣα原子の構築は１０回繰り返され、最小の目的関数値を持つＣα原子の座標が最適解として算出されるようにする。
ステップ６０：主鎖原子座標の構築及び最適化
ステップ５０（１）のＣαの原子座標に、他の主鎖原子（アミドＮ原子，カルボニルＣ原子，カルボニルＯ原子）及びＣαに化学結合したＣβ原子の座標を付加し、シミュレーティッドアニーリング法によって目的関数（後述の（８）式のＥ_ｍａｉｎ）を最小化するようにする。先ず、Ｃα原子の立体的な重ね合わせを行い、Ｃαの原子間距離が２．５Å以下の残基が取り上げられる。Ｃαを除く主鎖の原子座標は重なるべきＣα原子間距離が最小になるように参照タンパク質の座標から取得しモデル構造とする。
参照タンパク質の中に相当する残基が無い場合、主鎖の原子座標は、予め作成したデータベース（例えば、Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）中の相当する４残基のタンパク質断片から作成される。この過程の中で、残基ｉの主鎖原子はｉ−１番目からｉ＋２番目までのＣα原子間の最小のｒｍｓｄ値を持つ残基から選ばれる。その際、Ｎ末端の残基では、Ｃα原子座標の重ね合わせ範囲がｉ番目からｉ＋３番目までとなり、Ｃ末端の残基及びその一つ前の残基では同様にｉ−３番目からｉ番目まで及びｉ−２番目からｉ＋１番目までとなる。
主鎖原子の目的関数を基にシミュレーティッドアニーリング法によって主鎖原子座標（側鎖Ｃβ原子も含む。）が最適化される。目的関数は下記式（８）の通りである。従来法と本発明の重要な違いの一つは、当該式（８）中の化学結合ポテンシャルＥ_ｂｏｎｄと化学結合角ポテンシャルＥ_ａｎｇと化学結合ねじれ角ポテンシャルＥ_ｔｏｒ及び原子間相互作用ポテンシャルであるＥ_{ｎｏｎ−ｂｏｎｄ}を、ステップ４０で作成されたアライメント中の区切り文字Ｕを参考にして定めたｋ_Ｎの値（Ｎ番目のタンパク質鎖のＣ末端残基の通し番号に等しい。）によって下記のように修正した点である。

Ｅ_ｂｏｎｄは下記式（９）の通りである。

但し、
ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}かつｂ_ｉ＝Ｄ_{Ｃｉ，Ｎｉ＋１}の場合、加算しない。（ｉ＋１は添え字の添え字。）
ここで、ｂ_ｉ ^０は、標準の結合長であり、Ｎ−Ｃα、Ｃα−Ｃ、及びＣ−Ｎの３種類の化学結合によって異なるが、ここでは簡略化して表記している。Ｋ_ｂは定数であり、例えば２２５と設定される。Ｅ_ｂｏｎｄの計算は図５に示したように、ＣαにおけるＥ_ｌｅｎの計算と同様に、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間には化学結合が存在しないので、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}かつｂ_ｉ＝Ｄ_{Ｃｉ，Ｎｉ＋１}ならｂ_ｉ−ｂ_ｉ ^０＝０という条件により、Ｅ_ｂｏｎｄの計算に含めずに、相互作用の切断を行う。タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の区切り文字Ｕ_１からＵ_Ｍ−１までのＵの位置によって指定することができる。
Ｅ_ａｎｇは化学結合角の関数で、下記式（１０）の通りである。

但し、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}かつθ_ｉ＝θ_{Ｃαｉ，Ｃｉ，Ｎｉ＋１}の場合、加算しない。
ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}かつθ_ｉ＝θ_{Ｏｉ，Ｃｉ，Ｎｉ＋１}の場合、加算しない。
ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}かつθ_ｉ＝θ_{Ｃｉ，Ｎｉ＋１，Ｃαｉ＋１}の場合、加算しない。
（角度を決めている三つの原子の説明）
ここで、θ_ｉ ^０は、標準の結合角であり、それぞれの結合角の種類によって異なるが、ここでは簡略化して表記している。Ｋ_ａは定数であり、例えば４５と設定される。Ｅ_ａｎｇの計算は図６に示したように、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間には化学結合が存在しないので、Ｅ_ａｎｇの計算に含めずに、相互作用の切断を行う。タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の区切り文字Ｕ_１からＵ_Ｍ−１までのＵの位置によって指定することができる。
Ｅ_ｔｏｒは主鎖の化学結合ねじれ角ポテンシャルであり、下記式（１１）の通りである。

但し、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}ならφ_ｉ−φ_ｉ ^０＝ω_ｉ−ω_ｉ ^０＝φ_ｉ＋１−φ_ｉ＋１ ^０＝０
ここで、φ_ｉ ^０とφ_ｉ ^０は主鎖のねじれ角がＲａｍａｃｈａｎｄｒａｎプロットを満たすように設定する。即ち、（φ_ｉ，φ_ｉ）から最も座標が近く、かつＲａｍａｃｈａｎｄｒａｎプロットを満たすような（φ_ｉ ^０，φ_ｉ ^０）を選択する。また、ω_ｉ ^０は０としてｃｉｓ−Ｐｒｏ残基の場合のみπ（ｒａｄ）とする。Ｋ_ｔ及びＫ_ωは定数であり、例えば、それぞれ１０及び５０と設定される。Ｅ_ｔｏｒの計算は、Ｅ_ａｎｇの計算と同様に、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間には化学結合が存在しないので、Ｅ_ｔｏｒの計算に含めずに、相互作用の切断を行う。図７はωｉを示すものであるが、φ_ｉとφ_ｉ＋１の扱いも同様である。タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の区切り文字Ｕ_１からＵ_Ｍ−１までのＵの位置によって指定することができる。
Ｅ_{ｎｏｎ−ｂｏｎｄ}は原子間相互作用ポテンシャルで、下記式（１２）の通りである。

但し、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}のときｒ_ｉｊ∈｛ｒ_{Ｎｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｎｉ＋１}；ｒ_{Ｃｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｃβｉ＋１}；ｒ_{Ｃｉ，Ｃｉ＋１}；ｒ_{Ｏｉ，Ｎｉ＋１}；ｒ_{Ｏｉ，Ｃαｉ＋１}｝に対してｒ_ｉｊ≦８．０ならばＥ_{ｎｏｎ−ｂｏｎｄ}の計算に含める。
ここで、ε_ｉｊ ^＊とｒ_ｉｊ ^＊は定数で原子の種類によって異なる定数である（緒方浩二、東京理科大学学位論文，１９９９参照。）。Ｋ_ｎｏｎは定数であり、例えば０．２５と設定し、ｒ_ｉｊが、例えば８Å以下の原子対ｉ，ｊについて通常、それが３ボンド以上離れたものである場合には、Ｅ_{ｎｏｎ−ｂｏｎｄ}の計算に含めるものとする。Ｅ_{ｎｏｎ−ｂｏｎｄ}の計算は図８に示したように、ＣαにおけるＥ_ｖｄＷの計算と同様に、Ｎ番目のタンパク質鎖のＣ末端残基ｋ_ＮとＮ＋１番目のタンパク質鎖のＮ末端残基ｋ_Ｎ＋１の間には化学結合が存在しないので、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}のときは｛ｒ_{Ｎｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｎｉ＋１}；ｒ_{Ｃｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｃβｉ＋１}；ｒ_{Ｃｉ，Ｃｉ＋１}；ｒ_{Ｏｉ，Ｎｉ＋１}；ｒ_{Ｏｉ，Ｃαｉ＋１}｝に属するｒ_ｉｊが８Å以下であれば、新たにＥ_{ｎｏｎ−ｂｏｎｄ}の計算に含めなければならない。タンパク質鎖の総数がＭである場合には、この操作はＮ＝１からＭ−１まで行い、ｋ_１からｋ_Ｍ−１までの残基番号は、ステップ４０で作成されたアライメント中の区切り文字Ｕ_１からＵ_Ｍ−１までのＵの位置によって指定することができる。
Ｅ_ＳＳはＣｙｓ残基が生成するジスルフィド結合の関数で、下記式（１３）の通りである。

ここで、Ｄ_ｉ ^Ｃα及びＤ_ｉ ^Ｃβは、タンパク質鎖内及びタンパク質鎖間でジスルフィド結合を形成するＣｙｓ残基のペア同士のＣα間及びＣβ間距離である。残基番号ｉには全タンパク質鎖にわたる通し番号を用いるため、本発明ではタンパク質鎖間のＥ_ＳＳを従来法のポテンシャル関数で扱うことができる。Ｋ^ＳＳ _Ｃα及びＫ^ＳＳ _Ｃβは定数であり、例えば７．５である。
Ｅ_ｐｏｓは主鎖原子の位置に関する関数で、下記式（１４）の通りである。式の説明は前記式（６）に準ずる。

ここで〈ｗ_ｉｘ_ｉ〉は下記式（１５）のように与えられる。式の説明は前記式（７）に準ずる。

前記式（１２）の〈ｗ_ｉｘ_ｉ〉は、目的タンパク質及び参照タンパク質の間の構造の重ね合わせから求める。Ｋ_ｐｏｓは定数であり、例えば０．３である。
Ｅ_ｃｈｉはＣαの光学活性（キラリティー）に関するものであり、下記式（１６）の通りである。ここで、Ｃαのキラリティーとは、アミノ酸残基の光学異性体（Ｌ−体又はＤ−体）に関するものであり、通常はＬ−体のＣα原子になるように、当該式（１６）のポテンシャルが用いられる。

ここで、τ_ｉはｉ番目の残基のＮ−Ｃα−Ｃβ−Ｃで定められるねじれ角でありＫ_ｃｈｉは、例えば５０と設定される。
Ｅ_ｈｙｄｒは配列が類似したタンパク質中で保存された主鎖の水素結合に関するもので、下記式（１７）の通り定められる。

水素結合は、Ｎ原子とＯ原子の距離が２．９±０．５Åにあるときに設定される。複数の参照タンパク質中で水素結合があるか否かを判定するときは、７５％以上（四つ中三つ以上）の参照蛋白が存在すると認めた場合に水素結合ありと判定する。Ｋ_ｈｙｄｒは定数であり、例えば０．６である。
次にＣβを含む主鎖原子の最適化がシミュレーテッドアニーリングによって行われる。この過程で主鎖とＣβの原子の摂動が初期の位置に対して１．０Å以内になるようにする。この段階は主鎖とＣβの原子に対して通常２００回行われる。温度に相当するパラメータは通常５０又は２５から始まり一回毎に０．５倍にしてゆき０．０１になるまで続け、その後一定値とする。
主鎖の立体配置を幅広くサンプリングするために、好ましくは上記の方法を６回行い、最小の目的関数値（Ｅ_ｍａｉｎ）を持つ主鎖の原子座標を最適解とする。そして、温度に相当するパラメータは、通常最初の２回は５０からスタートして３回目から２５からスタートすることとする。
ステップ７０：側鎖原子座標の構築及び最適化
側鎖の構築は、大きく２段階に分かれており、「構造保存部位の側鎖構築」（ステップ７０（１）と「全体の側鎖構築」（ステップ７０（２））に分けられる（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。
ステップ７０（１）：構造保存部位の側鎖構築
算出された主鎖原子に対して、構造保存部位（ＳＣＲ領域）の場合、以前の研究における方法を用いて配列が類似したタンパク質から側鎖のねじれ角を得る。この方法の詳細は、″Ｔｈｅ ｒｏｌｅ ｏｆ ｐｌａｙｅｄ ｂｙ ｅｎｖｉｒｏｎｍｅｎｔａｌ ｒｅｓｉｄｕｅｓ ｉｎ ｓｉｄｅ−ｃｈａｉｎ ｔｏｒｓｉｏｎａｌ ａｎｇｌｅｓ ｗｉｔｈｉｎ ｈｏｍｏｌｏｇｏｕｓ ｆａｍｉｌｉｅｓ ｏｆ ｐｒｏｔｅｉｎｓ：Ａ ｎｅｗ ｍｅｔｈｏｄ ｏｆ ｓｉｄｅ ｃｈａｉｎ ｍｏｄｅｌｉｎｇ．″Ｏｇａｔａ Ｋ ａｎｄ Ｕｍｅｙａｍａ Ｈ，Ｐｒｏｔ．Ｓｔｒｕｃｔ．Ｆｕｎｃｔ．Ｇｅｎｅｔ．１９９８，３１，２５５−３６９に記載されている。この全内容は参考により本明細書中にその一部として組み込まれる。
この方法の中で配列が類似したタンパク質の中で保存されている側鎖の割合を算出し、この情報を基にして側鎖のモデリングを行う。側鎖の保存された部位の側鎖の原子座標は固定した主鎖原子に対して置かれる。例えば、配列が類似したタンパク質中でアルギニン残基のχ１角が保存されていれば、Ｃ_γ原子の座標を置くことができ、Ｐｈｅ残基でχ^１とχ^２角が保存されていれば、全ての側鎖原子を置くことができる。式（８）を用いたシミュレーティッドアニーリングの最適化の過程は、主鎖とＣβの原子のみ行われて、原子の摂動は１．０Å以内となるようにした。この主鎖とＣβの原子のアニーリングの段階は２００回行われる。そして、温度に相当するパラメータは２５からスタートして１回毎に０．５倍にしてゆき０．０１になるまで小さくなるようにする。前記式（８）の中のＥ_{ｎｏｎ−ｂｏｎｄ}は主鎖原子と部分的に作成された側鎖原子について行われる。そのとき側鎖原子の座標は最適化の過程を通じて保存されるようにする。
構造の情報である前記式（１４）のＭ_ｉと前記式（１７）の水素結合のＮ−Ｏのペアは最適化の過程で距離を計算し直し、特にＮ−Ｏペアについては、水素結合の有無判断を変化して用いられる。主鎖原子の配置を得るために、上記プロセスを３回繰り返し、目的関数の最小の主鎖原子の座標を算出構造とする。
ステップ７０（２）：全体の側鎖の構築
側鎖全体の構築は固定した主鎖及びＣβ原子の基で行う。これは上記したＯｇａｔａ Ｋ ａｎｄ Ｕｍｅｙａｍａ Ｈ，Ｐｒｏｔ．Ｓｔｒｕｃｔ．Ｆｕｎｃｔ．Ｇｅｎｅｔ．１９９８，３１，２５５−３６９に記載されている方法をもって行われ、それを用いることにより短時間で正確なモデルを与えることができる。先ず、前記式（８）の目的関数Ｅ_ｍａｉｎを用いて、主鎖構造（Ｃβを含む）を低温におけるモンテカルロ法によって最適化する。このとき、温度は０．００１に設定し、前記式（８）の中のＥ_{ｎｏｎ−ｂｏｎｄ}に関しては主鎖原子と側鎖全原子の間で計算を行う。そして、Ｎ、Ｃα、Ｃ、Ｃβ原子の最適化の過程で側鎖のねじれ角を最適化された状態を保つように側鎖の座標を再配置する。原子の摂動は０．５Å以内とする。次に、側鎖は削除され、上記の側鎖構築が繰り返される。このプロセスは２．４Åの原子同士のぶつかり合いがなくなり、かつＮ−Ｃα−Ｃβ−Ｃのねじれ角が−１２０±１５°の範囲に収まるまで繰り返される。
ステップ８０：モデル構造の妥当性の検査
でき上がったモデルの立体構造上の妥当性の検査は、ループの挿入や欠損部位を除いた領域で、モデルと参照タンパク質の主鎖原子座標の重ね合わせを行い、ｒｍｓｄが１Å以下であれば、妥当なモデルとした。
ステップ９０：最終構造の構築：立体構造予測
上記の通りステップ４０において得られたアライメントを基に、ステップ５０−８０において、今回、本発明者等が新に開発したｎｅｗＦＡＭＳ等のモデリングソフトを用いて立体構造を構築し、モデルが完成する。併せて、上記ステップ４０−８０で示す方法を「ｎｅｗＦＡＭＳ」と称する。一方、本発明の基本をなす従来型のモデリングソフト（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）を単に「ＦＡＭＳ」と称する。
尚、複数鎖モデリングの擬似単数鎖モデリング実現のために，単数鎖モデリングの式に対して前記式（２）、（３）、（４）、（９）、（１０）及び（１２）における本文中に示すような修正を加えたが、もう一つの方法として、それ等の式に対して該当する係数、前記式（２）のＫｌ、式（３）のＫａ、式（４）のＫｖｄｗ、式（９）のＫｂ、式（１０）のＫａ、式（１１）のＫｔを残基番号に依存する値として扱い、それぞれ調整することでも擬似単数鎖を実現できる。この擬似単数鎖モデリングは本発明を、具体的に実施する段階で、前記式（２）、（３）、（４）、（９）、（１０）及び（１２）の右辺の構成から、容易に考えられる。また、この擬似単数鎖モデリングは、後述の実施例１〜４と全く同じ結果を与える。
本発明によれば、配列上分離若しくは大きく間が開いた領域（複数の鎖）同士が空間的に接した立体構造を参照タンパク質に選択することもできる（図４４参照。）。配列上分離若しくは大きく間が開いた領域（複数の鎖）同士が空間的に接した状態の複合体の立体構造を、実験的に決定して（図４４ａ参照。）、これを参照タンパク質として本発明によりモデリングを行うと、目的タンパク質の立体構造を正確に構築することができる（本発明方法；図４４ｃ参照。）。これに対して、従来法により同様にモデリングを行うと、図４４ｂに示した通り両者の接触面が不正確なものが構築される。このような立体構造の構築においても、本発明方法によれば、従来法に比較してより正確に立体構造を構築することができることが理解される。
前記したように、本明細書中に引用される従来法についての文献の全内容は参考により本明細書中にその一部として組み込まれる。同様に、本願の基礎となる２００２年１月９日出願の日本出願：特願２００２−２８５９に含まれる明細書等出願書類の発明に関する全内容についても参考により本明細書中にその一部として組み込まれる。
好適な実施の形態
以下、実施例により本発明を詳細に説明する。これは、本発明を説明するためのものであり、故に本発明はこれ等実施例に限定されるものではない。
（実施例１） 代謝型グルタミン酸レセプター同族体のモデリング例
ラット由来の代謝型グルタミン酸受容体タイプＩ蛋白質（Ｍｅｔａｂｏｔｒｏｐｉｃ ｇｌｕｔａｍａｔｅ ｒｅｃｅｐｔｏｒ １；参照タンパク質）のアミノ酸１次配列は、Ｓｗｉｓｓ−Ｐｒｏｔより取得した。（エントリ名 ＭＧＲ１＿ＲＡＴ，アクセッション番号Ｐ２３３８５）。これは図１３に示す通り（配列表配列番号１参照。）、全体が１１９９残基で、前半の下線で示した４７８残基中、４４８〜４５６番目の９残基同士が接触部位となり、単量体（モノマー）が二つ集まった２量体（ホモダイマー）となり、グルタミン酸受容体部分を形成していることが、森川等によるＸ線構造解析の結果、明らかとなっている（Ｋｕｎｉｓｈｉｍａ，Ｎ．，Ｓｈｉｍａｄａ，Ｙ．，Ｔｓｕｊｉ，Ｙ．，Ｓａｔｏ，Ｔ．，Ｙａｍａｍｏｔｏ，Ｍ．，Ｋｕｍａｓａｋａ，Ｔ．，Ｎａｋａｎｉｓｈｉ，Ｓ．，Ｊｉｎｇａｍｉ，Ｈ．，Ｍｏｒｉｋａｗａ，Ｋ．：Ｓｔｒｕｃｔｕｒａｌ Ｂａｓｉｓ ｏｆ Ｇｌｕｔａｍａｔｅ Ｒｅｃｏｇｎｉｔｉｏｎ ｂｙ ａ Ｄｉｍｅｒｉｃ Ｍｅｔａｂｏｔｒｏｐｉｃ Ｇｌｕｔａｍａｔｅ Ｒｅｃｅｐｔｏｒ，Ｎａｔｕｒｅ ４０７，９７１，２０００参照。）。ＰＤＢには３種類のタンパク質立体構造（１ＥＷＫ，１ＥＷＴ，１ＥＷＶ）として登録されている。１ＥＷＫはリガンドとしてのグルタミン酸を含む構造であり、１ＥＷＴはリガンドとしてのグルタミン酸を含まない構造である。Ｘ線構造解析の結果からは、１ＥＷＫはモノマー同士が比較的接近したｃｌｏｓｅｄ形であり、一方、１ＥＷＴはモノマー同士が比較的離れたｏｐｅｎ形であることが分かっている。両者は、このようにモノマー同士の相対配置が異なるために、ドメイン間の接合面の状態も大きく異なっている。
即ち、グルタミン酸受容体タイプＩ蛋白質では、リガンドが結合していない状態（１ＥＷＴ）ではｏｐｅｎ形であるが、リガンドが結合した状態（１ＥＷＫ）ではｃｌｏｓｅｄ形が安定化することが推察されている。従って、このタンパク質では、モノマー同士の相対配置を含む接合面の状態を正確にモデリングすることが機能解明のために非常に重要である。そこで、実施例としては、類似配列を探索するのに１ＥＷＴを参照タンパク質として、リガンドを含まないタンパク質２量体の立体構造モデル（参照タンパク質）を構築した。また、従来法でダイマーを構成する個々のモノマーをモデリングし、両者を併せてダイマーモデルとした場合と本発明による方法でダイマー自体をモデリングした場合で、ダイマー間の接合面の状態に、立体構造上の優劣（エネルギー的な安定性の差）が生じるかを比較した。
先ず、ＭＧＲ１＿ＲＡＴのアミノ酸配列１１９９残基をクエリーとして、２００１年１１月時点のＰＩＲをモチーフプロファイルとし、ＰＩＲ，Ｓｗｉｓｓ Ｐｒｏｔ，ＴＲＥＭＢＬ，ＴＲＥＭＢＬ＿ＮＥＷ，ＧｅｎＰｅｐｔ等のタンパク質アミノ酸配列データベース（何れも２００１年１１月時点）の７７４８０４配列に対してＰＳＩ−ＢＬＡＳＴ検索を行った。Ｅ−ｖａｌｕｅが０．００１以下の条件で実施した結果、１４５０９個の類似配列とのアライメントを取得した。そのうち、Ｅ−ｖａｌｕｅが０となるものが７０ある（ホモロジーは２３〜１００％で、Ｅ−ｖａｌｕｅが非常に小さいため計算機出力としてゼロと記述されている。）。これ等は、機能的にほぼ同一と見てよい。この中で同じＲＡＴ由来のレセプターでＳｗｉｓｓ Ｐｒｏｔ中のエントリー名「ＭＧＲ５＿ＲＡＴ」アクセッション番号「Ｐ３１４２４」のタンパク質２量体（目的タンパク質）のモデリングを行った（Ａｂｅ Ｔ．，Ｓｕｇｉｈａｒａ Ｈ．，Ｎａｗａ Ｈ．，Ｓｈｉｇｅｍｏｔｏ Ｒ．，Ｍｉｚｕｎｏ Ｎ．，Ｎａｋａｎｉｓｈｉ Ｓ．，Ｍｏｌｅｃｕｌａｒ ｃｈａｒａｃｔｅｒｉｚａｔｉｏｎ ｏｆ ａ ｎｏｖｅｌ ｍｅｔａｂｏｔｒｏｐｉｃ ｇｌｕｔａｍａｔｅ ｒｅｃｅｐｔｏｒ ｍＧｌｕＲ５ ｃｏｕｐｌｅｄ ｔｏ ｉｎｏｓｉｔｏｌ ｐｈｏｓｐｈａｔｅ／Ｃａ２＋ｓｉｇｎａｌ ｔｒａｎｓｄｕｃｔｉｏｎ，Ｊ．Ｂｉｏｌ．Ｃｈｅｍ．２６７，１３３６１−１３３６８，１９９２参照。）。ＭＧＲ５＿ＲＡＴはイノシトールリン酸とカルシウムイオンが関係する情報伝達タンパク質で、ラット由来の代謝型グルタミン酸受容体サブタイプＶタンパク質（ｍｅｔａｂｏｔｒｏｐｉｃ ｇｌｕｔａｍａｔｅ ｒｅｃｅｐｔｏｒ ｓｕｂｔｙｐｅ５）であり、アミノ酸残基数は１２０３個である。参照配列ＭＧＲ１＿ＲＡＴと目的配列ＭＧＲ５＿ＲＡＴの間のホモロジーは６２．２％であり、図１４に両者のアライメント結果を示す（配列表配列番号１及び２参照。）。
（ＭＧＲ１＿ＲＡＴとＭＧＲ５＿ＲＡＴのアライメント結果）

図１４のアライメント中でＰＤＢの１ＥＷＴの座標が存在する部分のみを表示すると、図１５のようになる（配列表配列番号３及び４参照。）。このアライメントを用いてモデリングを行った。
（ＭＧＲ１＿ＲＡＴとＭＧＲ５＿ＲＡＴのアライメント結果；１ＥＷＴで座標が存在する部分のみ）

図１５のアライメントはモノマーのものであり、従来型ＦＡＭＳ（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）の入力ファイル形式では、図１６のようになる（配列表配列番号３及び４参照。）。このモデリングは従来型ＦＡＭＳを使用することによって実行した。モデリング結果を図１７に示した。
更に、ホモダイマーとして扱った場合は、図１６のアライメントを「Ｕ」の文字でつなぐことになり、図１８の通りである。この「Ｕ」を使う入力ファイル形式は、本発明者等によって開発されたものである。モデリング結果を図１９に示した。
また、図１９のｎｅｗＦＡＭＳによるホモダイマーモデルの接合面の拡大図を図２０に示す。この図において、タンパク質−タンパク質相互作用の認識部位が主鎖や側鎖で２．４Å以内の衝突はなく精度よくモデリングされている。一方、１ＥＷＴのホモダイマーの各々のモノマーの３次元座標を使って、独立に従来型ＦＡＭＳでＭＧＲ５＿ＲＡＴの各々の立体構造を構築すると、接合面では２．４Å以内の原子間接触が８ヶ所で生じた。この構造の場合、接合面は衝突が多いことによりエネルギー的には非常に不安定である。一方、ｎｅｗＦＡＭＳによるモデル構造では、接合面で衝突がないことによりエネルギー的に安定である。このことは本発明者等が開発した複数鎖モデリングを行うｎｅｗＦＡＭＳの優れていること（新規性）を示すものである。
（実施例２） トランスポーターのモデリング例
多剤耐性の原因の一つと考えられているＡＢＣトランスポーターの同族体の立体構造がＸ線結晶回折で４．５Åの分解能で解析されておりＣαの座標のみではあるがＰＤＢ：１ＪＳＱで登録されている（参照タンパク質）。この立体構造によってＡＢＣトランスポーターは、細胞膜上でホモダイマーを形成し、モノマー間の相対位置のフリップフロップ運動からリン脂質を細胞から放出する働きを持つことが示唆されている（Ｇｅｏｆｆｒｅｙ Ｖｈａｎｇ ａｎｄ Ｃｈｉｓｔｏｐｈｅｒ Ｂ．Ｒｏｔｈ，ＳＣＩＥＮＣＥ，Ｖｏｌ ２９３，ｐｐ．１７９３参照。）。
１ＪＳＱは８本鎖（Ａ鎖、Ｂ鎖、Ｃ鎖、Ｄ鎖、Ｅ鎖、Ｆ鎖、Ｇ鎖、Ｈ鎖）あり、Ａ−Ｂ鎖、Ｃ−Ｄ鎖、Ｅ−Ｆ鎖、Ｇ−Ｈ鎖のホモダイマー４組の形で登録されている。本実施例では、先ず８個の鎖をＣαの座標のみからそれぞれ従来型ＦＡＭＳの自動モデリングにより主鎖及び側鎖を構築した。８個の側鎖まで構築した座標を見て、Ｂ鎖の主鎖に関する化学結合ねじれ角φ及びφが最もラマチャンドランプロット上のエネルギー的に不安定な領域を侵していない構造であった。次に、キメラモデリング法（Ｔ．Ｙｏｎｅｄａ，Ｈ．Ｋｏｍｏｏｋａ，Ｈ．Ｕｍｅｙａｍａ，Ｊ．Ｐｒｏｔ．Ｃｈｅｍ．，１６，５９７−６０５，１９９７参照。）を用いてＢ鎖を基本構造として、部分的にその他の鎖で補いモノマーの構造をモデリングした。次に、Ｂ鎖のモノマーの座標を回転・移動しＡ鎖と重ね合わせ、移動後の座標を新たなＡ鎖とした。このようにして得たＡ−Ｂ鎖のペアを鋳型構造ＭＳＢＡとして、下記のモデリングに用いている（モデル化した参照タンパク質）。１ＪＳＱのモノマーのアミノ酸配列は５５５残基であり（２ケ所の下線部は座標欠損で、これ等を除くと４５０残基となる）、図２１の通りである（配列表配列番号５参照。）。
この５５５残基の配列を照会する目的配列（クエリー）として、２００１年１１月時点のＰＩＲデータベースを使用してクエリーのモチーフプロファイルを作成し、ＰＩＲ，Ｓｗｉｓｓ Ｐｒｏｔ，ＴＲＥＭＢＬ，ＴＲＥＭＢＬ＿ＮＥＷ，ＧｅｎＰｅｐｔ等のタンパク質アミノ酸配列データベース（何れも２００１年１１月時点）の７７４８０４配列に対してＩＭＰＡＬＡ検索（ＰＳＩ−ＢＬＡＳＴと似ているが、Ｓｍｉｔｈ ａｎｄ Ｗａｔｅｒｍａｎ法によるアライメントを用いる手法：Ａ．Ａ．Ｓｃｈａｆｆｅ ｅｔ ａｌ．，ＢＩＯＩＮＦＯＲＭＡＴＩＣＳ，１５（１２），１０００−１０１１，１９９９参照。）を行った。Ｅ−ｖａｌｕｅが０．００１以下の条件で実施した結果、１３７０５個の類似配列とのアライメントを取得した。
一例として、ＳｗｉｓｓＰｒｏｔデータベース中のＩＤ「Ｏ９３４３７」に対してモデリングを行った。Ｏ９３４３７（目的タンパク質）はＳｗｉｓｓ−Ｐｒｏｔのホームページにおいて、ニワトリＡＢＣトランスポータータンパク質と記述されている（Ｅｄｅｌｍａｎｎ Ｈ．Ｍ．Ｌ．，Ｄｕｃｈｅｋ Ｐ．，Ｒｏｓｅｎｔｈａｌ Ｆ．Ｅ．，Ｆｏｅｇｅｒ Ｎ．，Ｇｌａｃｋｉｎ Ｃ．，Ｋａｎｅ Ｓ．Ｅ．，Ｋｕｃｈｌｅｒ Ｋ．，”Ｃｍｄｒ１，ａ ｃｈｉｃｋｅｎ Ｐ−ｇｌｙｃｏｐｒｏｔｅｉｎ，ｃｏｎｆｅｒｓ ｍｕｌｔｉｄｒｕｇ ｒｅｓｉｓｔａｎｃｅ ａｎｄ ｉｎｔｅｒａｃｔｓ ｗｉｔｈ Ｅｓｔｒａｄｉｏｌ”，Ｂｉｏｌ．Ｃｈｅｍ．３８０，２３１−２４１，１９９９参照。）。アミノ酸残基数は１２８８個で、多剤耐性を示し、卵胞ホルモンであるエストラジオールと相互作用するタンパク質である。Ｏ９３４３７をモデリングするため、前述のＭＳＢＡとアライメントした結果を図２２に示した（配列表配列番号６及び７参照。）。
（「Ｏ９３４３７」とＭＳＢＡとのアライメント結果）

但し、参照タンパク質であるＭＳＢＡには、途中に大きな座標欠損が２箇所あるため（図２１の１ＪＳＱのアミノ酸配列中、下線で示した部分）、図２２のアライメントのままではモデリングができない。そこで、座標欠損の位置に「Ｕ」の文字を挿入し、アライメントを図２３の様に修正した。これは、一つのタンパク質をあたかも三つのタンパク質から成るように記述し、複数タンパク質のモデリングを実行するものである。
（「Ｏ９３４３７」とＭＳＢＡの新しいアライメント結果）

これをｎｅｗＦＡＭＳの入力である入力ファイル形式にすると、図２４に示したようになる。
図２４のアライメントはタンパク質を擬似単数鎖として扱って、タンパク質モデルを構築しているが、この入力ファイルはあくまでモノマーのものである。このアライメントを用いたモデリング結果を図２５に示した。
次に、ホモダイマーとしてモデル構築をする場合は、上記のアライメントを更に「Ｕ」の文字でつなぐことになり以下の図２６のようになる。
図２６の入力ファイル形式を用いたｎｅｗＦＡＭＳによるモデリング結果を図２７に示した。
上記ホモダイマーの接触部分のグラフィック表示は図２８の通りである。
このように座標欠損がある場合も擬似的に複数鎖として扱い「Ｕ」の文字で置き換えてやることにより、ファンデアワールスの原子の衝突を考慮した構造を算出することが可能となり、相互作用の認識部位を正確に表したホモダイマーをモデル構築できる。作成したホモダイマーモデル「Ｏ９３４３７」を観察すると、Ｘ線解析よりホモダイマーの膜貫通部位と考えられる部位で、トランスポートされる物質が結合の可能性があると述べられている場所は、図２９の丸で囲まれた領域であることが推察された。この領域は非常に多くの疎水性アミノ酸残基が存在し、「Ｏ９３４３７」の機能とされている卵胞ホルモンエストラジオールが容易に結合できることが推察される。因みに、相互作用するエストラジオールは図３０に示されるように、分子全体が疎水性基から構成されている。一方、上述の１ＪＳＱに記述されたＣα原子座標のみを基本にして主鎖と側鎖をモデル構築したホモダイマーの各々モノマーの３次元座標を参照タンパク質にして、独立に「Ｏ９３４３７」ＡＢＣトランスポーターを各々モデル構築して、相互作用の認識部位を観察すると、ホモダイマー間の相互作用が考慮されていないため、接合面では２．４Å以内の原子間接触が５４ヵ所で生じた。上記のようにリガンドのエストラジオールが結合する可能性のある近くに２量体接触領域があることから、このホモダイマーの接合部位のモデル構築は機能を正確に記述するために極めて大切である。このことはやはり本発明における複数鎖モデリングの方法の優れた点（新規性）や、有用性を示しているものと考えられる。
（実施例３） 炭疸菌の致死因子タンパク質に結合したペプチドの改変例
２００１年１１月号の英国雑誌ＮＡＴＵＲＥ（Ｐａｎｎｉｆｅｒ ｅｔ．ａｌ，ＮＡＴＵＲＥ，ｖｏｌ ４１４，ｐｐ．２２９−２３３）に近年世間を騒がせている炭疸菌の毒素である致死因子（ＬＦ、分子量９万）の立体構造がＸ線結晶回折法により３．９０Åの分解能で報告されている。このタンパク質は、四つのドメインから成り炭疸の病原性に不可欠な酵素であり、ヒト細胞の一つ又は複数の情報伝達系を阻害する。具体的には、ｍｉｔｇｅｎ−ａｃｔｉｖａｔｅｄ ｐｒｏｔｅｉｎ ｋｉｎａｓｅ ｋｉｎａｓｅ（ＭＡＰＫＫ）と呼ばれるタンパク質ファミリーのＮ末端部分と接触してその末端部分を切断すると報告されている。ＬＦは極めて特異性の高いプロテアーゼである。ＰＤＢでは、ＬＦ単体が１Ｊ７Ｎとして、ＭＡＰＫＫ−２と呼ばれるＭＡＰＫＫファミリーのＮ末端部１６残基とＬＦの複合体が１ＪＫＹ（参照タンパク質）として登録されている。ＭＡＰＫＫ−２のＮ末端１６残基はＬＦのドメインＩＩ，ＩＩＩ，ＩＶの三つから作られる長く深い溝にはさまって切断される。この毒素をターゲットにした薬は有望であると考えられているので、ＭＡＰＫＫ−２のＮ末端１６残基の立体構造を参考にペプチドではない薬物を設計することが可能である。
本実施例では、このＭＡＰＫＫ−２の１６残基のアミノ酸残基に仮想的に変異を起こさせた配列（目的タンパク質）のモデルを構築し、ＬＦとの相互作用の様子がどのように変わるかの一例を示す。通常、薬物の設計には薬物の腸管吸収等を考慮する必要から疎水性の結合部位に焦点を定める。ＰＤＢ：１ＪＫＹはＡ鎖（ＬＦ）とＢ鎖（ＭＡＰＫＫ−２の１６残基）があり、アミノ酸配列は図３１の通りである（配列表配列番号８及び９参照。）。
このＢ鎖のＭＡＰＫＫ−２の１６残基のうち、第二残基と第三残基「ＬＡ」の部分は、Ａ鎖のＬＦとの間に空間的に隙間が空いている（図３２参照）。そこで、この部分のアミノ酸残基をより嵩高い側鎖を持つアミノ酸残基に置き換え疎水性相互作用の接触をより強固にするような阻害剤を設計することにより、その阻害剤がＭＡＰＫＫ−２の１６残基のＮ末端と競合的に結合することによって、ＬＦのプロテアーゼとしての活性を阻害できる可能性がある。
「ＬＡ」の換わりに、例えばこの嵩高い側鎖をもつアミノ酸残基として「ＦＦ」を選んだ場合のモデル構築を行った。
このＡ鎖（ＬＦ）とＢ鎖（ＭＡＰＫＫ−２）の複合体（１ＪＫＹ）を鋳型として、図３３の入力ファイル形式を用いて（配列表配列番号１０参照。）、ｎｅｗＦＡＭＳによりモデル構築を行った。モデリング結果を図３４に示した。
以上のように「ＬＡ」を「ＦＦ」に変えてモデル構築した一例を示したが、ｎｅｗＦＡＭＳを用いてそれ以外の変異体についても網羅的に解析することが可能であり、その中からより効果的な阻害剤のリード化合物を推定することが可能になる。
（実施例４） 実際にＰＤＢに登録されている複数鎖から成るタンパク質立体構造と構築モデルの比較による本発明の検証例（ブラインド−テスト）
プロティナーゼの一種であるトリプシンとそのタンパク質阻害剤である膵臓トリプシンインヒビター（ＰＴＩ）との複合タンパク質２ＰＴＣの立体構造はＸ線結晶解析法により、１．９０Åの分解能で解かれて登録されている（Ｍａｒｑｕａｒｔ，Ｍ．，Ｗａｌｔｅｒ，Ｊ．，Ｄｅｉｓｅｎｈｏｆｅｒ，Ｊ．，Ｂｏｄｅ，Ｗ．，Ｈｕｂｅｒ，Ｒ．：Ｔｈｅ Ｇｅｏｍｅｔｒｙ ｏｆ ｔｈｅ Ｒｅａｃｔｉｖｅ Ｓｉｔｅ ａｎｄ ｏｆ ｔｈｅ Ｐｅｐｔｉｄｅ Ｇｒｏｕｐｓ ｉｎ Ｔｒｙｐｓｉｎ，Ｔｒｙｐｓｉｎｏｇｅｎ ａｎｄ ｉｔｓ Ｃｏｍｐｌｅｘｅｓ ｗｉｔｈ Ｉｎｈｉｂｉｔｏｒｓ，Ａｃｔａ Ｃｒｙｓｔａｌｌｏｇｒ．，Ｓｅｃｔ．Ｂ ３９，４８０，１９８３参照。）。２ＰＴＣはＥ鎖のトリプシンとＩ鎖のＰＴＩからなる複合体で登録されている。セリンプロテアーゼであるトリプシンは膵臓から不活性なトリプシノーゲンとして分泌され十二指腸液中の酵素Ｅｎｔｅｒｏｋｉｎａｓｅにより活性型トリプシンになる。トリプシンはＥｎｄｐｅｐｔｉｄａｓｅ一種でＡｒｇやＬｙｓのような塩基性アミノ酸のカルボキシル基側のペプチド結合を切断するタンパク分解酵素である。
一方、セリンプロテアーゼの一種であるトロンビンＥ１９２Ｑ酵素（グルタミン酸１９２がグルタミンに置換されている。）とＰＴＩの複合タンパク質１ＢＴＨの立体構造はＸ線結晶解析法により、２．３Åの分解能で解かれて登録されている（ｖａｎ ｄｅ Ｌｏｃｈｔ，Ａ．，Ｂｏｄｅ，Ｗ．，Ｈｕｂｅｒ，Ｒ．，Ｌｅ Ｂｏｎｎｉｅｃ，Ｂ．Ｆ．，Ｓｔｏｎｅ，Ｓ．Ｒ．，Ｅｓｍｏｎ，Ｃ．Ｔ．，Ｓｔｕｂｂｓ，Ｍ．Ｔ．，“Ｔｈｅ ｔｈｒｏｍｂｉｎ Ｅ１９２Ｑ−ＢＰＴＩ ｃｏｍｐｌｅｘ ｒｅｖｅａｌｓ ｇｒｏｓｓ ｓｔｒｕｃｔｕｒａｌ ｒｅａｒｒａｎｇｅｍｅｎｔｓ：ｉｍｐｌｉｃａｔｉｏｎｓ ｆｏｒ ｔｈｅ ｉｎｔｅｒａｃｔｉｏｎ ｗｉｔｈ ａｎｔｉｔｈｒｏｍｂｉｎ ａｎｄ ｔｈｒｏｍｂｏｍｏｄｕｌｉｎ”ＥＭＢＯ Ｊ．１６，２９７７，１９９７参照。）。このタンパク質はＬとＨ鎖からなるトロンビンＥ１９２ＱとＰ鎖のみから成るＰＴＩ、同じように、Ｊ鎖とＫ鎖から成るトロンビンＥ１９２ＱとＱ鎖のみから成るＰＴＩのホモダイマー２組の形で登録されている。ここでは、前者の１組のタンパク質複合体に注目する。トロンビンは血液凝固因子の一つであるプロトロンビンの活性型蛋白であり、このトロンビンが形成されることによってフィブリノーゲンがフィブリンとなり血液凝固を起こす。
本実施例では、本発明ソフトであるｎｅｗＦＡＭＳで２ＰＴＣのＥ鎖とＩ鎖から１ＢＴＨのＨ鎖とＰ鎖のモデリングを行い、実際にＰＤＢ（Ｐｒｏｔｅｉｎ Ｄａｔａ Ｂａｎｋ）に登録されている１ＢＴＨのＸ線構造との比較を行った。但し、１ＢＴＨモデルの挿入ループの構築に使用するＣα原子や主鎖原子のデータベースの中には、１ＢＴＨ自身は含まれておらず、モデルの精度を検証するに当たって、正解構造を含むデータベースは使用していないことを確認してある。また、従来型のＦＡＭＳでの１ＢＴＨ＿Ｈ、１ＢＴＨ＿Ｅを個々にモデリングしたものとの比較もｎｅｗＦＡＭＳの精度を検証するために示す。ここで、注目すべきは、１ＢＴＨのＨ鎖はトロンビンの１９２番目のアミノ酸残基がグルタミン酸からグルタミンに残基置換されていることである。トロンビンは１９２番目のアミノ酸残基はグルタミン酸であり、通常、酵素のループ部分が立体的に阻害してＰＴＩとは結合しない。しかし、１９２番目のアミノ酸残基をグルタミンに置換する（トロンビンＥ１９２Ｑ）とＰＴＩと結合するときに、上記の妨害ループ部分がＰＴＩの影響を受けて適合誘導で動き、両者のタンパク質は結合する。それが１ＢＴＨの名で登録された複合タンパク質である。
目的タンパク質である１ＢＴＨのＨ鎖及びＰ鎖のアミノ酸配列は図３５に示す通りである（配列表配列番号１１及び１２参照。）。１ＢＴＨのＨ鎖のアミノ酸残基数は２５４残基であり、Ｐ鎖のアミノ酸残基数は５８残基である。
一方、参照タンパク質である２ＰＴＣのＥ鎖及びＩ鎖のアミノ酸配列は図３６に示す通りである（配列表配列番号１２及び１３参照。）。２ＰＴＣのＥ鎖のアミノ酸残基数は２２３残基であり、Ｉ鎖のアミノ酸残基数は５８残基である。
１ＢＴＨのＨ鎖と２ＰＴＣのＥ鎖におけるアライメントは図３７に示す通りである。ホモロジーは４３．５％であった。
（１ＢＴＨのＨ鎖と２ＰＴＣのＥ鎖におけるアライメント結果）

１ＢＴＨのＰ鎖と２ＰＴＣのＩ鎖におけるアライメントは図３８に示す通りである。ホモロジーは１００％であった。
（１ＢＴＨのＰ鎖と２ＰＴＣのＩ鎖におけるアライメント結果）

各アライメントの結果の従来型のＦＡＭＳの入力ファイル形式は図３９の通りである。
今回のｎｅｗＦＡＭＳの入力ファイル形式では、上記のアライメントを「Ｕ」の文字でつなぐことにより、図４０の通りである。
今回ｎｅｗＦＡＭＳによってモデル構築をＨ鎖とＰ鎖に限って行った結果、実際にＰＤＢに登録されている１ＢＴＨと今回のｎｅｗＦＡＭＳによる１ＢＴＨモデルとのｒ．ｍ．ｓ．ｄ値はＨ鎖、Ｐ鎖共に含め全体で２．１１Åとなった。また、実際にＰＤＢに登録されている１ＢＴＨと従来のＦＡＭＳによる１ＢＴＨモデルとの全体のｒ．ｍ．ｓ．ｄ値は２．１８Åであった。ｒ．ｍ．ｓ．ｄ値は若干ながら改善されていることが分かる。
次に、従来のＦＡＭＳによるモデルのアミノ酸の原子同士が衝突しているかどうか、Ｈ鎖とＰ鎖の接触面で２．４Å未満のアミノ酸残基接触があるかを検証した。その結果、Ｈ鎖とＰ鎖の分子間では７ヵ所で原子間接触が生じていた。
一方、今回の本発明で使用可能なソフトｎｅｗＦＡＭＳによるモデルのＨ鎖とＰ鎖間でアミノ酸の原子同士が衝突しているかどうかを観察すると、２．４Å以下は無かった。
以上で示したように、従来型のＦＡＭＳでは、原子同士で衝突を起こしてしまうが、今回のｎｅｗＦＡＭＳでは、原子同士では衝突は発見できなかった。セリンプロテアーゼである１ＢＴＨのＨ鎖はＰ鎖との接触領域において、活性部位である触媒部位と基質結合部位をもち、機能の正確な記述のためには、これ等Ｈ鎖とＰ鎖の接触領域を精度良くモデル構築せねばならない。この実施例はＸ線解析があり、それと比較しながら、複数鎖の擬似単数鎖モデル構築が精度の高いことを示した。このことは、重ねてｎｅｗＦＡＭＳの優れていること（新規性）を示したことになる。
（実施例５）
閲覧したいタンパク質を、生物種固有の任意の記号、タンパク質コード名、参照タンパク質名、１行程度の機能説明の文字列との部分一致等の、論理積検索によって、簡便に事前の知識無く目的とするタンパク質にたどり着けるよう設計されたインターフェースを開発することができた。図４３に、一例として、ＧＴＯＰ（国立遺伝学研究所で公開しているゲノムのアミノ酸配列と立体構造データベースＰＤＢのアミノ酸配列とのアライメントのホームページ）で閲覧できるアライメントに基いて構築した立体構造モデルデータベースへのインターフェース画面を示す。これは、３次元座標に基いた閲覧したいタンパク質立体構造が上記の単数鎖若しくは複数鎖のタンパク質立体構造データベース中に存在する場合、生物種固有の任意の記号、タンパク質コード名、参照タンパク質名、１行程度の機能説明の文字列との部分一致等の、論理積検索によって、簡便に事前の知識無く目的とするタンパク質にたどり着けるよう設計されたインターフェースとなっている。
この例の場合、４１生物種のａｅｒｏ，ａｆｕｌ，，，，，等の記号はＧＴＯＰで命名されたものに基く（２００１年９月時点）。立体構造、即ち３次元座標に基いた閲覧したいタンパク質がある場合、その生物種の記号の横のチェックボックスにチェックを入れておく。更に、タンパク質コード名、参照タンパク質名、１行程度の機能説明等の論理積によって絞込み検索機能を持っている。
発明の効果
本発明によれば、ペプチドを含むタンパク質について、複数鎖から成る立体構造未知のタンパク質（これをコードする遺伝子を含む。）の立体構造予測又は複数鎖の各々の成分のアミノ酸残基が改変された複合タンパク質の立体構造予測に関して、従来法より効率的かつ単純に、信頼性が高い立体構造情報を得ることができる。その結果、ゲノム配列解析、ＤＮＡチップを用いた発現プロファイリング解析、プロテオーム解析等で興味深い遺伝子やタンパク質が見つかった場合に、そのタンパク質の立体構造を基準とした機能予測を行うことが可能となる。これに伴ってそのペプチドやタンパク質の効率的な改変が可能となる。また、タンパク質の機能部位を予測することにより、従来法に比べて、より効率的又は効果的にタンパク質や低分子化合物の新規医薬品を設計するための情報を取得することができる。
産業上の利用分野
以上、本発明の方法は、医農薬の分子設計等を中心に、生物情報の解析を行う分野（バイオインフォマティックス）において極めて有用であると考えられる。単数鎖のアミノ酸配列に適応できる蛋白質自動モデリングシステムの適用範囲を複数鎖から成るタンパク質に拡張すると同時に、リガンド分子、受容体分子或いは酵素等の多様な改変を可能にするものであり、当該発明の有用性の更なる増大が期待される。
本発明は、産業上多くの分野、特に医薬品、食品、化粧品、医療、構造解析等の分野で広く実施することができ、故に極めて有用である。
【配列表】

【図面の簡単な説明】
［図１］ 図１は、本発明の説明のために使用するタンパク質の全アミノ酸配列と各鎖Ｃ末端残基の通し番号ｋ_Ｎの関係を示したものである。
全鎖数＝Ｍ、各ｋ_Ｎは区切り文字Ｕより識別する。Ｎはタンパク質鎖の番号を表す。
［図２］ 図２は、ポテンシャル計算式のＥｌｅｎ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである（Ｃαの場合のみ）。
［図３］ 図３は、ポテンシャル計算式のＥａｎｇ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである（Ｃαの場合のみ）。
［図４］ 図４は、ポテンシャル計算式のＥｖｄｗ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである（Ｃαの場合のみ）。
［図５］ 図５は、ポテンシャル計算式のＥｂｏｎｄ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである。
［図６］ 図６は、ポテンシャル計算式のＥａｎｇ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである。
［図７］ 図７は、ポテンシャル計算式のＥｔｏｒ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである。
［図８］ 図８は、ポテンシャル計算式のＥｎｏｎ−ｂｏｎｄ項におけるＣ末端及びＮ末端の取り扱い方を図示したものである。
［図９］ 図９は、本発明による複数鎖から成るタンパク質の立体構造構築法の一例を示すフローチャートである。
［図１０］ 図１０は、目的タンパク質のＣα原子座標の構築方法の一例を図示したものである。
アライメントの一致部分は参照タンパク質から取得し、無い部分はＮ、Ｃ両端それぞれ重なった２残基の重ね合わせのｒｍｓｄが最小のものをデータベースから取得する（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。
［図１１］ 図１１は、ローカルスペースホモロジー（ＬＳＨ）の概念を図示したものである。
例えば、図の円（球）の中心のＴ残基に関する計算では、灰色の残基で考慮される。アライメントにおける四角で囲った部分が考慮される残基ペアであり、＊のマークがあるところとの比率がＬＳＨ（この例では５６．２％）である（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。
［図１２］ 図１２は、ローカルスペースホモロジー（ＬＳＨ）と構造保存部位（ＳＣＲｓ）にある比率との関係を図示したものである。
ＬＳＨは目的タンパク質と参照タンパク質とのＣα原子の重ね合わせから計算され、ＳＣＲｓにある比率は目的タンパク質の全残基数に対するＳＣＲｓ中の残基数を表す（Ｋ．Ｏｇａｔａ ａｎｄ Ｈ．Ｕｍｅｙａｍａ，Ｊｏｕｒｎａｌ ｏｆ Ｍｏｌｅｃｕｌａｒ Ｇｒａｐｈｉｃｓ ａｎｄ Ｍｏｄｅｌｉｎｇ １８，２５８−２７２，２０００参照。）。
［図１３］ 図１３は、実施例１で使用するラット由来の代謝型グルタミン酸受容体タイプＩの全アミノ酸配列を示したものである。
［図１４］ 図１４は、図１４−１と図１４−２で構成され、実施例１において、ＭＧＲ１＿ＲＡＴとＭＧＲ５＿ＲＡＴのアライメント結果を示したものである。
［図１５］ 図１５は、実施例１において、ＭＧＲ１＿ＲＡＴとＭＧＲ５＿ＲＡＴのアライメント結果を示したものである。
１ＥＷＴで座標が存在する部分のみ示されている。
［図１６］ 図１６は、実施例１において、従来型ＦＡＭＳの入力ファイル形式を示したものである。
［図１７］ 図１７は、実施例１において、ｍＧｌｕレセプター Ｓｗｉｓｓ Ｐｒｏｔ「ＭＧＲ５＿ＲＡＴ」のモノマーモデルを示したものである。
［図１８］ 図１８は、実施例１において、ｎｅｗＦＡＭＳに対する入力ファイル形式を示したものである。
［図１９］ 図１９は、実施例１において、ｍＧｌｕレセプターＳｗｉｓｓ Ｐｒｏｔ「ＭＧＲ５＿ＲＡＴ」のダイマーモデルを示したものである。
［図２０］ 図２０は、実施例１において、ｍＧｌｕレセプターＳｗｉｓｓ Ｐｒｏｔ「ＭＧＲ５＿ＲＡＴ」ダイマーモデルの接合面の拡大図である。
［図２１］ 図２１は、実施例２において、１ＪＳＱのモノマーのアミノ酸配列を示したものである。
［図２２］ 図２２は、実施例２において、「Ｏ９３４３７」とＭＳＢＡとのアライメント結果を示したものである。
［図２３］ 図２３は、実施例２において、「Ｏ９３４３７」とＭＳＢＡの新しいアライメント結果を示したものである。
［図２４］ 図２４は、実施例２において、「Ｏ９３４３７」に対するｎｅｗＦＡＭＳの入力ファイル形式を示したものである。
［図２５］ 図２５は、実施例２において、ＡＢＣトランスポーターＳｗｉｓｓＰｒｏｔ「Ｏ９３４３７」のモノマーモデルを示したものである。
［図２６］ 図２６は、実施例２において、ホモダイマーに対してｎｅｗＦＡＭＳで用いる入力ファイル形式を示したものである。
［図２７］ 図２７は、実施例２において、ＡＢＣトランスポーターＳｗｉｓｓＰｒｏｔ「Ｏ９３４３７」のダイマーモデルを示したものである。
［図２８］ 図２８は、実施例２において、ＡＢＣトランスポーターＳｗｉｓｓＰｒｏｔ「Ｏ９３４３７」のダイマーモデル接合面を図示したものである。
［図２９］ 図２９は、実施例２において、ＡＢＣトランスポーターモデルにおける疎水性領域を図示したものである。
［図３０］ 図３０は、実施例２において、エストラジオールの分子構造を示したものである。
［図３１］ 図３１は、実施例３において、ＰＤＢ：１ＪＫＹのアミノ酸配列を示したものである。
ＰＤＢ：１ＪＫＹはＡ鎖（ＬＦ）とＢ鎖（ＭＡＰＫＫ−の１６残基）を有する。
［図３２］ 図３２は、実施例３において、ＬＦ中のＭＡＰＫＫ−２を図示している。
尚、Ｌｅｕ−２とＡｌａ−３についてはボールアンドスティックモデルを示す。
［図３３］ 図３３は、実施例３において、「ＬＡ」から「ＦＦ」へのペプチド改変におけるｎｅｗＦＡＭＳの入力ファイル形式を示したものである。
［図３４］ 図３４は、実施例３において、ＬＦ中のＭＡＰＫＫ−２改変体を図示している。
尚、Ｐｈｅ−２とＰｈｅ−３についてはボールアンドスティックモデルを示す。
［図３５］ 図３５は、実施例４において、１ＢＴＨのＨ鎖及びＰ鎖のアミノ酸配列を示したものである。
［図３６］ 図３６は、実施例４において、２ＰＴＣのＥ鎖及びＩ鎖のアミノ酸配列を示したものである。
［図３７］ 図３７は、実施例４において、１ＢＴＨのＨ鎖と２ＰＴＣのＥ鎖におけるアライメントを図示したものである。
［図３８］ 図３８は、実施例４において、１ＢＴＨのＰ鎖と２ＰＴＣのＩ鎖におけるアライメントを図示したものである。
［図３９］ 図３９は、実施例４において、従来型のＦＡＭＳの入力ファイル形式を示したものである。
［図４０］ 図４０は、実施例４において、ｎｅｗＦＡＭＳの入力ファイル形式を示したものである。
［図４１］ 図４１は、実施例４において、１ＢＴＨのＨ鎖とＰ鎖のモデル接合面を図示したものである。
［図４２］ 図４２は、実施例４において、２ＰＴＣからモデリングされた１ＢＴＨ（Ｈ鎖とＰ鎖）モデルと実際にＰＤＢに登録されている１ＢＴＨ（Ｈ鎖とＰ鎖）の立体構造の重ね合わせを示したものである。
濃い方が実際のＸ線結晶解析法によるものであり、薄い方が今回モデリングを行ったものである。
［図４３］ 図４３は、実施例５において、立体構造データベースへのインターフェース画面の一例を図示したものである。
３次元座標閲覧可能なモデルの一覧が表示され、右のボタンをクリックすることによりアライメントが表示され、３次元座標を取得できるように構成されている。
［図４４］ 図４４は、配列上分離若しくは大きく間が開いた領域（複数の鎖）同士が空間的に接した立体構造を参照タンパク質に選択して、本発明により立体構造をモデリングした場合と、従来法によりモデリングした場合とを概念的に図示して比較したものである。
図４４ａ：配列上分離若しくは大きく間が開いた領域（複数の鎖）同士が空間的に接した状態の複合体の立体構造を、実験的に決定したもの；図４４ｂ：図４４ａに示した構造を参照タンパク質として従来法によりモデリングしたもの（好ましくないモデル）；図４４ｃ：図４４ａに示した構造を参照タンパク質として本発明によりモデリングしたもの（好ましい良いモデル）。

Claims (7)

1. 与えられた任意のアミノ酸配列を有する単数鎖より成るタンパク質の立体構造を構築するホモロジーモデリング法又はスレッディング法を拡張して（拡張型モデリング法）、与えられた任意のアミノ酸配列を有する複数鎖より成る目的タンパク質の立体構造を、コンピュータが以下のステップを実行することにより構築する方法であって、
   （ａ）複数鎖を有する目的タンパク質のアミノ酸配列を用意し、用意された各タンパク質鎖のアミノ酸配列に対して、コンピュータが参照タンパク質の立体構造データベースから参照タンパク質を検索し、目的タンパク質と検索した複数の参照タンパク質とのアミノ酸配列アライメントを行うアライメントステップ、
   （ｂ）前記アライメントの結果に基づいて、コンピュータが目的タンパク質の各タンパク質鎖のアミノ酸配列と一致度の高い一つ又は複数の参照タンパク質のアミノ酸配列を選択するステップ、
   （ｃ）目的タンパク質及び選択された参照タンパク質のそれぞれに含まれる各タンパク質鎖のアミノ酸配列のＣ末端部に区切り文字を付加して擬似単数鎖としての形態を有するように修正する配列アライメントの修正ステップ、及び
   （ｄ）前記修正されたアライメントの情報に基づいて、目的タンパク質主鎖アミノ酸中のＣα原子について選択された参照タンパク質の立体構造から座標を取得し、目的関数を最小化するようにＣαの原子座標を最適化し、最適化されたＣαの原子座標に主鎖の他の原子座標（Ｃβの原子座標を含む。）を付加して目的関数を最小化するように主鎖の原子座標を最適化し、最適化された主鎖の原子座標に側鎖の他の原子座標を付加し目的関数を最小化するように側鎖の原子座標を最適化するステップ
   を含み、
   前記工程（ｄ）において目的タンパク質の立体構造を構築するための目的関数としてのポテンシャル計算式は、前記修正後の擬似単数鎖から各タンパク質鎖の末端残基番号を判定し、Ｃα原子座標の構築においてはその境界におけるＣα原子間距離ポテンシャル（擬似化学結合ポテンシャル）及びＣα原子結合角ポテンシャル（擬似化学結合角ポテンシャル）の切断、並びに当該境界におけるＣα原子間ファンデルワールスポテンシャル（原子間相互作用ポテンシャル）の付加を行うことで、主鎖原子座標及び側鎖原子座標の構築においてはその境界における化学結合ポテンシャル、化学結合角ポテンシャル及び化学結合ねじれ角ポテンシャルの切断、並びに当該境界における原子間相互作用ポテンシャルの付加を行うことで、選択された参照タンパク質の立体構造データに基づきＣα原子座標、主鎖原子座標及び側鎖原子座標の構築を行い、シミュレーティッドアニーリング法、分子力学計算、及びモンテカルロ法の少なくとも１種により擬似エネルギー値を表す目的関数の最小化（最適化）を行う
   ことを特徴とする複数鎖を有するタンパク質の立体構造構築方法。
2. 前記工程（ａ）において、参照タンパク質を検索しアライメントを出力するソフトウエアがＦＡＭＳ、ＦＡＳＴＡ、ＰＳＩ−ＢＬＡＳＴ、ＬＩＢＲＡ、ＲＢＳ−ＢＬＡＳＴ、ＩＭＰＡＬＡ、ＣｌｕｓｔａｌＷ、ＨＭＭＥＲ、及びＢＩＯＣＥＳの少なくとも１種である請求項１記載の方法。
3. 目的タンパク質を構成する複数鎖のうち少なくとも二つの鎖がポリペプチド鎖であり、タンパク質−タンパク質相互の認識部位に位置するアミノ酸残基のそれぞれについて２０種類のアミノ酸の可能な組合せに基いて、ポテンシャルエネルギー値を指標として、適合性の優れた、又は劣った改変体のデータセットを作成して、当該少なくとも二つのポリペプチド鎖の各々のタンパク質の機能を上昇させた、又は低下させた立体構造を構築できるようにした請求項１記載の方法。
4. 目的タンパク質を構成する複数鎖のうち少なくとも一つの鎖がアミノ酸誘導体又はペプチド誘導体（ペプチド性リガンド）であり、かつ参照タンパク質において該当するリガンド分子と類似した化学構造を持つものである場合に、目的タンパク質の当該誘導体を新たな残基名及び１文字コードで定義し、かつ参照タンパク質の当該リガンドを更に新たな残基名及び１文字コードで定義したアライメントを手動的又は自動的に作成して、当該リガンドの配列を構成するそれぞれの残基について２０種類のアミノ酸及びその誘導体の可能な組み合わせに基づいて、ポテンシャルエネルギー値が低い順にランク付けした場合の上位複数を、受容体タンパク質の結合部位への適合性の優れた改変体として格納したアミノ酸誘導体又はペプチド誘導体リガンドモデルデータセットを構築できるようにした請求項１記載の方法。
5. 目的タンパク質を構成する複数鎖の成分の少なくとも一つがペプチド性リガンドであり、当該リガンドのアミノ酸配列を固定し、当該リガンドを認識する部位に位置するアミノ酸残基のそれぞれについて２０種類のアミノ酸の可能な組み合わせに基づき、ポテンシャルエネルギー値を指標として、上位複数の受容体タンパク質の結合部位への適合性の優れた改変体のデータセットを作成して、当該リガンドに結合し得る多様なリガンド受容体タンパク質の立体構造を構築できるようにした請求項１記載の方法。
6. 前記工程（ｂ）において、目的タンパク質の各タンパク質鎖のアミノ酸配列に対してＥ−ｖａｌｕの小さい参照タンパク質のアミノ酸配列を選択する請求項１記載の方法。
7. 前記目的タンパク質の立体構造の構築におけるポテンシャル計算式が下記の内容を含む請求項１記載の方法：
   全鎖数＝Ｍのときのポテンシャル計算式については、Ｎがタンパク質鎖の番号を、ｋ_ＮがＮ番目のタンパク質鎖におけるＣ末端残基の通し番号をそれぞれ表し、ｉ＝１，．．．，Ｍ−１を簡略化して、ｉ＝１，Ｍ−１と表記するとき、
   （Ａ）Ｃα原子座標の構築と最適化プロセスにおける計算では、擬似化学結合ポテンシャルのｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の場合を含めず、擬似化学結合角ポテンシャルのｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}＋１の場合を含めない。また、原子間相互作用ポテンシャルの場合は、ｉ＝ｋ_Ｎ−１ならｊ＞ｉ＋１、ｉ＝ｋ_Ｎならｊ＞ｉを各々追加する。
   （Ｂ）主鎖原子座標の構築と最適化プロセスにおける計算では、化学結合ポテンシャルにはｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}のときのＣｉとＮｉ＋１間の結合を含めず、化学結合角ポテンシャルには、Ｃ、Ｏをそれぞれカルボニルの炭素原子と酸素原子、Ｃαをα炭素原子、Ｎを窒素原子とするとき、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の角度Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１、Ｏｉ−Ｃｉ−Ｎ_ｉ＋１及びＣｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１を含めず、更に、化学結合ねじれ角ポテンシャルには、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}のときの角度Ｎｉ−Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１、Ｃαｉ−Ｃｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１、及びＣｉ−Ｎ_ｉ＋１−Ｃα_ｉ＋１−Ｃ_ｉ＋１を含めない。また、原子間相互作用ポテンシャルには、原子間距離をｒとするとき、ｉ＝ｋ_{Ｎ（Ｎ＝１，Ｍ−１）}の、ｒ_ｉｊ∈｛ｒ_{Ｎｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｎｉ＋１}；ｒ_{Ｃαｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｎｉ＋１}；ｒ_{Ｃｉ，Ｃαｉ＋１}；ｒ_{Ｃｉ，Ｃβｉ＋１}；ｒ_{Ｃｉ，Ｃｉ＋１}；ｒ_{Ｏｉ，Ｎｉ＋１}；ｒ_{Ｏｉ，Ｃαｉ＋１}｝についてｒ_ｉｊ≦規定値の場合を追加する。

Images