説明

分子間の類似度を評価するための高速グラフマッチ検索装置及び方法

【課題】分子グラフに関して、2分子間の原子対応を求め該対応に基づいて2分子を重ね合わせする方法を実現する高速グラフマッチアルゴリズムを提供する。
【解決手段】第1の分子Aの夫々の原子と第2の分子Bの夫々の原子を対応付け(m(Ai)=Bk)て重ね合わせ、ABの間の最適な原子間対応を求める方法であって、S1(Ai,Bk)、S2(Ai,Bk)及びS3(Ai,Bk)の各々を全ての{i,k}の組について求める工程と、その工程にて最大のS3(Ai,Bk)を算出した原子(Ai,Bk)の対応から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の組が無くなるまで続けたときの全体の対応におけるグラフマッチスコアM(A,B)を求める工程と、その工程におけるM(A、B)が閾値より大きいならば、AとBにつき重ね合わせを出力する工程を含む方法を提供する。

【発明の詳細な説明】
【技術分野】
【0001】
本発明は、高速グラフマッチ検索アルゴリズムを利用して、2分子間の原子対応を求め対応に基づいて2分子を仮想的に重ね合わせ、2分子間の類似度を求めて評価する、高速グラフマッチ検索装置及び方法に関する。
【背景技術】
【0002】
医薬や農薬の分子設計において、2つの分子に係る分子構造を仮想空間にて重ね合わせすることが頻繁に行われる。図13は、そのような、2つの分子(Cholic acid[CHD]とCorticosteron[COR])を仮想空間にて重ね合わせすることを模式的に示す図である。しかしながら、2つの分子についての最適な重ね合わせを探索し決定することは非常に困難な問題である。
【0003】
例えば、分子Aと分子Bとの重ね合わせの問題について、片方の分子Aを『CMP』とした場合に、それに基づき重ね合わせにて探索可能な重ね合わせの対象の分子Bを求める場合を検討する。ここで「探索可能な」というのは、全探査を8時間労働・週休2日の労働時間で50年程度行って、探索が解決され得ると想定される、という程の意味である。例えば、人手による計算による場合では、分子BがCysteinである場合、1.3×10通り程度の重ね合わせの計算を行い、最適な重ね合わせを求めることが可能となる(図14(a))。同様に、デスクトップコンピュータによる場合では、分子BがDiaminopimelateである場合、1.5×1015通り程度の重ね合わせの計算を行い、最適な重ね合わせを求めることが可能となる(図14(b))。更に同様に、超高速度電子計算機による場合でも、分子BがAMPである場合、8.3×1021通りの重ね合わせの計算を行い、最適な重ね合わせを求めることが可能となる(図14(c))。このように、分子Aと分子Bの最適な重ね合わせを全探索に拠ることは、膨大な時間が掛かるため、必ずしも現実的な方法ではない。
【0004】
よって、2分子間の原子対応を求め該対応に基づいて2分子の最適な重ね合わせを実現するグラフマッチにおいて、多少の間違いを許容しつつも発見的に高速に行うことが求められている。
【先行技術文献】
【特許文献】
【0005】
なお、化合物検索のアルゴリズムに関する先行技術文献として、以下のような6件が挙げられる。
【0006】
【特許文献1】特許第4001657号
【特許文献2】特許第3928000号
【特許文献3】国際出願01/097094号
【特許文献4】国際出願02/41184号
【特許文献5】国際出願2007/004643号
【非特許文献】
【0007】
【非特許文献1】J.Computer-Aided Molecular Design, 13:499-512, 1999 Estimation of active confirmations of drugs by a new molecular superposing procedure
【発明の概要】
【発明が解決しようとする課題】
【0008】
本発明は、原子をノード、化学結合をエッジとして表現した分子グラフに関して、2分子間の原子対応を求め該対応に基づいて2分子を重ね合わせする方法を高速に実現する、グラフマッチ検索装置及び方法を提供することを目的とする。
【課題を解決するための手段】
【0009】
本発明は、上記の目的を達成するために為されたものである。本発明に係る請求項1に記載の、分子間の類似度を評価するための高速グラフマッチ検索装置は、
第1の分子Aを構成する原子(Ai,Aj,・・・)の各々に係る座標データと第2の分子Bを構成する原子(Bk,Bl,・・・)の各々に係る座標データを記憶部から入力し、演算部にロードされるコンピュータプログラムに従って、演算部及び記憶部に構築される仮想メモリ空間において第1の分子Aの夫々の原子(Ai,Aj,・・・)と第2の分子Bの夫々の原子(Bk,Bl,・・・)との対応付け(m(Ai)=Bk)を求めて重ね合わせを行い(i,j,k,lはいずれも自然数)、第1の分子Aと第2の分子Bの間の最適な原子間対応、及び第1の分子Aと第2の分子Bの類似度に係るデータを出力部に出力する、第1の分子Aと第2の分子Bとの類似度を評価するための高速グラフマッチ検索装置において、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、周囲の環境が相互にどれだけ似ているかを示す第1の類似指標S1(Ai、Bk)を求める第1の算出手段と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、等しい結合距離にある周囲の原子Aj、Blの全ての組につき、第1の類似指標S1(Aj,Bl)を積算して算出する第2の類似指標S2(Ai、Bk)を求める算出手段であって、その原子Ai、Bkの対の各原子から等しい結合距離にある周囲の原子Aj、Blが同じ元素であれば、更に第1の類似指標S1(Aj,Bl)に係数を掛けた上で積算する、第2の類似指標S2(Ai、Bk)を求める第2の算出手段と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対を始点とし、第1の分子Aの原子と第2の分子Bの原子とを順次対応付けして全体の対応を作成し、そのときに算出されるグラフマッチスコアM(A,B)を値とする第3の類似指標S3(Ai、Bk)を求める算出手段であって、対応付け作成時には、既に対応付け済みの原子に直接結合する原子を次に選択すること、及び第2の類似指標S2が高い対を選択することを優先することを、条件とする、第3の類似指標S3(Ai、Bk)を求める第3の算出手段と、
第3の算出手段にて最大のS3(Ai,Bk)を算出した際の、始点の原子(Ai,Bk)の対から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の組が無くなるまで続けたときの、全体の対応におけるグラフマッチスコアM(A,B)を求める第4の算出手段と、
第4の算出手段におけるグラフマッチスコアM(A、B)が閾値より大きいならば、第1の分子Aと第2の分子Bにつき第4の算出手段で算出した原子間対応及びグラフマッチスコアM(A,B)を出力する第5の出力手段と
を含むことを特徴とする。
【発明の効果】
【0010】
本発明により、原子をノード、化学結合をエッジとして表現した分子グラフに関して、2分子間の原子を対応させ該対応に基づいて2分子を重ね合わせするにあたり、最適な重ね合わせを高速に且つ精度よく求めることができる。
【図面の簡単な説明】
【0011】
【図1】本発明の実施形態に係るグラフマッチによる分子構造の高速アルゴリズムを実現するコンピュータシステムの構成の例を示す図である。
【図2】本発明の実施形態に係る高速グラフマッチ探索アルゴリズムによる分子構造の重ね合わせ及びその表示のためのプログラムのフローチャートである。
【図3】分子Aと分子Bの、原子(ノード)及び結合(エッジ)を模式的に示す図(図3(1))と、図3(1)に示す分子Aと分子Bに基づいて算出された分子グラフマッチスコアの例(図3(2))である。
【図4】図2に示すステップS10において、分子Aと分子Bの間の原子対応関係{m(Ai)}とグラフマッチスコアM(A、B)を求める高速グラフマッチ探索アルゴリズムのフローチャートである。
【図5】分子Aの一部及び分子Bの一部を示す図であって、原子Aiと原子Bkに関する指標S1の算出を説明するための図である。
【図6】分子Aの一部及び分子Bの一部を示す図であって、原子Aiと原子Bkに関する指標S2の算出を説明するための図である。
【図7】分子Aの一部及び分子Bの一部を示す図であって、原子Aiと原子Bkに関する指標S3の算出を説明するための図である。
【図8】分子Aの一部及び分子Bの一部を示す図であって、図4のステップS1008におけるグラフマッチスコアM0の算出を説明するための図である。
【図9】分子Aの一部及び分子Bの一部を示す図であって、図4のステップS1012における微調整を説明するための図である。
【図10】ねじれ角を調整して、分子Aと分子Bをより重ね合わせて表示することを模式的に示す図である。
【図11】クエリの原子{Ai}に対応した{m(Ai)}の組で、重ね合わせを行い原子座標を出力した図(図11(a))と、クエリ構造から共通骨格にあたる原子座標を出力した図(図11(b))である。
【図12】総組み合わせ数1012以下の問題に対して、全探査により最大スコアを求め、本発明の実施形態に係る高速グラフマッチ探索アルゴリズムによる解と比較を行った際の、全探査組み合わせ数に対する計算時間をグラフ化したもの(図12(1))と、全探査組み合わせ数に対する正解率をグラフ化したもの(図12(2))と、正解スコア差と累積正解率の関係をグラフ化したもの(図12(3))である。
【図13】2つの分子(Cholic acid[CHD]とCorticosteron[COR])を仮想空間にて重ね合わせすることを模式的に示す図である。
【図14】分子Aと分子Bとの重ね合わせの問題について、片方の分子Aを『CMP』とした場合に、それに基づき重ね合わせにて探索可能な重ね合わせの対象の分子Bの例を示した図である。
【発明を実施するための形態】
【0012】
以下、図面を参照して本発明に係る好適な実施の形態を説明する。
【0013】
本実施形態に係るグラフマッチによる分子構造の高速アルゴリズムは、コンピュータを用いて行われるものであり、C言語などの適切なプログラム言語によって記述されたプログラムをコンピュータで実行し、(後で説明する)様々な分子を構成する原子に関する座標データをコンピュータ上で構築される仮想メモリ空間に展開することにより、実現されるものである。
【0014】
図1は、本実施形態に係るグラフマッチによる分子構造の高速アルゴリズムを実現するコンピュータシステム2の構成の例を示す図である。コンピュータシステム2は、ディスプレイ等の出力部12、キーボード16やマウス18などの入力部、並びに、演算部、記憶部及び通信制御部等を含む中央処理部14から構成される。中央処理部14は、インターネット4等の外部ネットワークを介して、外部サーバ8や外部データベース10と接続しそれら外部サーバ8や外部データベース10とデータを送受信することができるように、構成されている。
【0015】
本実施形態で利用される、様々な分子についての原子座標に係るデータは、PDB(プロテインデータバンク;蛋白質構造データバンク)フォーマットのデータであり、通常、外部の商用及び公開データベース10等から提供される。例えば、PDBフォーマットの様々な分子についての原子座標に係るデータは、外部ネットワーク4を介して外部の商用及び公開データベース10からダウンロードされ、コンピュータシステム2に付属する記憶部に格納される。これらのデータは、図2に示すフローチャートに係る処理を実行する際、記憶部から読み出されて利用される。
【0016】
1.高速グラフマッチ探索アルゴリズムによる分子構造の重ね合わせ処理
図2は、本実施形態に係る高速グラフマッチ探索アルゴリズムによる分子構造の重ね合わせ及びその処理のフローチャートである。図2を参照して本実施形態に係る分子構造の重ね合わせ処理を説明する。まず、重ね合わせの一方の分子(分子Aとする)についてのPDBフォーマットの原子座標を読み込む(ステップS02)。読み込んだPDBフォーマットの原子座標に基づいて、分子Aの結合距離・結合次数・回転可能結合の設定を行う(ステップS04)。分子の結合距離・結合次数・回転可能結合の設定については後で説明する。
【0017】
ステップS02及びステップS04と並行して、重ね合わせのもう一方の分子(分子Bとする)についてのPDBフォーマットの原子座標を読み込む(ステップS06)。なお、分子Bは複数であることがある。次に、分子Bの一つについて結合距離・結合次数・回転可能結合の設定を行う(ステップS08)。
【0018】
続いて、高速グラフマッチ探索アルゴリズムを行い、分子Aの原子(Ai)から分子B(Bk)への対応関係{m(Ai)}及びそのときのグラフマッチスコアM(A、B)を求める(i、kはいずれも自然数)(ステップS10)。ここで、グラフマッチスコアM(A、B)とは、2分子間の原子対応を求め該対応に基づいて2分子の最適な重ね合わせを実現するグラフマッチにおいて、最適さの程度を示す指標である。なお、グラフマッチスコアM(A、B)、対応関係{m(Ai)}、及び高速グラフマッチ探索アルゴリズムの、夫々の詳細については、後で説明する。
【0019】
グラフマッチスコアM(A、B)が閾値より大きいかどうか確認される(ステップS12)。閾値より大きいということは、そのグラフマッチスコアM(A、B)を実現する重ね合わせのための対応関係(m(Ai))が十分に適切であることを意味する(ステップS12のYes)。このとき、分子Aに対する分子Bのねじれ角が調節され(ステップS14)、分子Aと分子Bにつき原子アラインメント及び構造重ね合わせが出力される(ステップS16)。ねじれ角の調節、並びに、原子アラインメント及び構造重ね合わせの出力についても、後述する。
【0020】
更に、次の分子Bがあるかどうか判断される(ステップS18)。次の分子Bがあれば(ステップS18のYes)、次の分子Bについての結合距離・結合次数・回転可能結合の設定(ステップS08)以降の処理が繰り返される。
【0021】
分子Bが無くなれば(ステップS18・No)、出力部12に基本骨格(又は共通骨格)を出力して(ステップS20)処理を終了する。
【0022】
2.結合距離・結合次数・回転可能結合の設定
図2のステップS04及びS08で行われる「結合距離・結合次数・回転可能結合の設定」について説明する。
【0023】
(2.1)結合距離
PDBフォーマットに係るデータが示す分子構造では、原子間の結合が定義されていないことがある。そこで本実施形態では、一つの分子において、原子iと原子jの間の原子間距離が2.00Åより短い場合は化学結合が存在するものとしてデータ上、化学結合を設定する(i、jはいずれも自然数)。この「原子間距離」は、PDBから読み込まれる原子座標に基づいて計算される。更に、一つの分子において二つの原子を取り上げたとき、それら2原子を繋ぐ化学結合の数を「結合距離」とする。それら2原子を繋ぐ経路が複数存在するときは最小のものを取る。結合を一つずつ延長することで、一つの分子内の全ての原子間に結合距離が設定される。
【0024】
(2.2)結合次数
PDBフォーマットに係るデータが示す分子構造では、原子間の結合次数が定義されておらず、且つ、一般に水素原子を含んでいない。そこで、以下の表1の示すルールに従い、原子間距離に基づき結合次数を求める。
【表1】

【0025】
(2.3)回転可能結合
直接結合する原子の対(原子iと原子j)の全てについて、上記「(2.1)結合距離」の定義プロセスを、原子の対間の結合が存在しないものとして実行する。その結果、原子の対(原子iと原子j)間に結合距離が設定されず、且つ、原子iと原子jの間の結合が単結合である場合は、原子iと原子jの対の間の結合は「回転可能結合」であると設定する。
【0026】
3.分子グラフマッチスコア定義
本実施形態に係るグラフマッチによる分子構造の高速アルゴリズムでは、分子グラフマッチスコアM(A,B)を定義している。なお{M(A,B)}は、分子Aと分子Bとの間の分子グラフマッチスコアであることを示す。図3(1)は、分子Aと分子Bの、原子(ノード)及び結合(エッジ)を模式的に示す図である。
以下に、本実施形態で利用する分子グラフマッチスコアM(A,B)の定義((定義1)、(定義2)、(定義3)及び(定義4))について説明する。
【0027】
(定義1);「Ai」は、分子Aのi番目の原子であることを示す。「Ai−Aj」は、AiとAjの結合を示す。
【0028】
(定義2);分子Aの原子i(Ai)が、分子Bの原子k(Bk)に対応することを、「 m(Ai)=Bk 」と表すものとする。即ち、m(Ai)=Bkとは、分子Aの原子iが対応する分子Bの原子kを示す。
【0029】
(定義3)
分子グラフマッチスコアM(A,B)は以下の式(数1)で定義される
【数1】

数1の各項は、以下の通り定義される。なおE(Ai,Aj)は、実行modeにより異なる値を持つ。この実行modeは、図1に示す入力部等を介して事後的に外部から設定され得るものである。
【表2】

【0030】
図3(2)は、上述の定義に従い、図3(1)に示す分子Aと分子Bに基づいて算出された分子グラフマッチスコアの例である。模様が同じであれば同じ元素であり、エッジは全て単結合であるとしているので、実行modeに関わり無く、図3(2)に示す値(特に、M(A、B)=14)となる。
【0031】
4.高速グラフマッチ探索アルゴリズム
図4は、図2に示すステップS10において、分子Aと分子Bの間の、原子の対応関係{m(Ai)}とグラフマッチスコアM(A、B)を求める高速グラフマッチ探索アルゴリズムのフローチャートである。以下、このフローチャートを参照し、高速グラフマッチ探索アルゴリズムを具体的に説明する。
【0032】
[ステップS1002];まず、分子Aを構成する原子と、分子Bを構成する原子との全ての組み合わせ(Ai,Bk)について、以下の数2及び表3で定義される「S1(Ai,Bk)」を求める。
【数2】

【表3】

【0033】
S1(Ai,Bk)は、原子Aiと原子Bkの対において、周囲の環境(同じ結合距離に同じ種類の原子があるか)がどれだけ似ているかを示す指標である。
【0034】
例えば、図5に示される分子Aの一部、及び分子Bの一部において、原子Aiから2の結合距離にある原子Ajと、原子Aiに対応する原子Bkから2の結合距離にある原子Blとが同一元素であれば、s1(Aj,Bl)の値は“1”になる。原子Aiと原子Bkの対を中心として、同じ結合距離にある、分子Aの原子と分子Bの原子が同じかどうか、全{j,l}の組について確認し、“1”又は“0”を設定して積算する。上記のS1(Ai,Bk)は、対応する原子Ai,Bkからみて、同じ結合距離の位置に同じ元素がある、という場合が多い程、大きくなる。
【0035】
[ステップS1004];次に、以下の数3及び表4で定義される「S2(Ai,Bk)」を、全ての{i,k}の組について求める。
【数3】

【表4】

【0036】
S2(Ai,Bk)は、対応する原子Ai、Bkの対の夫々において、その対の各原子から等しい結合距離にある周囲の原子Aj、Blの全ての組について、上記の、周囲の環境がどれだけ似ているかを示す指標であるS1(Aj,Bl)を積算する指標であるが、その対の各原子から等しい結合距離にある周囲の原子Aj、Blが同じものであれば、更にS1(Aj,Bl)に係数(上記表では12)を掛けて積算される。従って、対応する原子Ai、Bkの対の各原子について、周囲の環境が類似し、更に周囲の環境のその周囲の環境が類似すれば、大きくなる指標である。
【0037】
例えば、図6に示される、原子Aiを含む分子Aの一部、及び原子Bkを含む分子Bの一部において、S2(Ai,Bk)を検討する。原子Aiからある結合距離(図6では2)にある原子Ajと、原子Bkからそれと等しい結合距離にある原子Blとの全ての対につき、s2(Aj,Bl)、即ちS1(Aj,Bl)、又はS1(Aj,Bl)×12を積算する。特に、AjとBlが同じ元素であれば、S1(Aj,Bl)は所定数倍(ここでは12倍)されて積算されて、S2が求められる。原子Ai及び原子Bkからの結合距離は、1から最大値(即ち、原子Ai又は原子Bkから最も遠い原子までの結合距離)まで変動することが想定される。上述のとおり、S1(Aj,Bl)は、原子Aj,Blの対において、(図6のAm、Bnなどの)周囲の環境がどれだけ似ているかを示す指標である。
【0038】
上記のS2(Ai,Bk)では、対応する2つの原子Ai,Bkに関して、同じ結合距離の位置の原子の対(Aj,Bl)のS1(Aj,Bl)が積算されるが、(Aj,Bl)が同じ元素であれば、S1(Aj,Bl)が所定数倍(12倍)されて積算されるから、周囲の原子の構成が近似するように対応付けされていると、やはりS2(Ai,Bk)は大きくなる。なお、係数「12」は別の数値であってもよい。
【0039】
[ステップS1006];次に、以下の数4及び表5で定義される「S3(Ai,Ak)」を、全ての{i,k}の組について求める。
【数4】

【表5】

【0040】
S3(Ai,Bk)は、原子Ai、Bkの対を始点とし、次々に分子Aの原子と分子Bの原子を対応付けして全体の対応を作成し、そのときのグラフマッチスコアM(A,B)を値とする指標である。ここで、対応付け作成時には、既に対応付け済みの原子に直接結合する原子を次に選択すること、及び指標S2が高い対を選択するのを優先することを、条件としている。
【0041】
例えば、図7に示される、原子Aiを含む分子Aの一部、及び原子Bkを含む分子Bの一部において、S3(Ai,Bk)を検討する。始点は、原子Ai、Bkの対である。原子Aiには、原子Aj、Ap、Arが直接結合する。原子Bkには、原子Bl、Bq、Bsが直接結合する。{Aj、Ap、Ar}と{Bl、Bq、Bs}とから形成され得る原子同士の(3×3=9通りの)対のうちから、(Aj、Bl)の対のS2が最大であるとすると、原子Ajと原子Blを対応付けすることになる。
【0042】
次に、分子Aにおいて対応付けが済んだAi−Ajには、原子Ap、Ar、At、Avが直接結合する。分子Bにおいて対応付けが済んだBk−Blには、原子Bq、Bs、Bu、Bwが直接結合する。{Ap、Ar、At、Av}と{Bq、Bs、Bu、Bw}とから形成され得る原子同士の(4×4=16通りの)対のうちから、(Ap、Bq)の対のS2が最大であるとすると、原子Apと原子Bqを対応付けすることになる。これにより、分子Aにおいては、原子Ai、Aj、Apの対応付けが完了し、分子Bにおいては、原子Bk、Bl、Bqの対応付けが完了する。
【0043】
このような対応付けを、対応可能原子の対が無くなるまで、順次繰り返して行う。対応付けが終われば、その対応付けの下でのグラフマッチスコアMを求める。このような対応付け及びグラフマッチスコアM算出が、全ての{i,k}の組について行われる。
【0044】
上記のS3(Ai,Bk)では、全ての{Ai,Bk}の組み合わせの各々において、原子の対の始点{Ai,Bk}の周囲から徐々に、S2(対応する原子Ai、Bkの対について、周囲の環境が類似し、更に周囲の環境のその周囲の環境が類似すれば、大きくなる指標)の大きさに着目して、分子Aの原子と分子Bの原子とが対応付けされ、グラフマッチスコアが計算されることになる。
【0045】
[ステップS1008];次に、ステップS1006にて最大のS3(Ai,Bk)を算出した際の、始点の原子(Ai,Bk)の対応から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の対が無くなるまで続け、全体の対応におけるグラフマッチスコアM0(A,B)を求める。このとき、途中、原子の対応の対と、次の原子の対応の対とにおいて、分子Aの原子は直接結合していなくてもよく、同様に、分子Bの原子も直接結合していなくてもよい。
【0046】
例えば、図8に示される、原子Aiを含む分子Aの一部、及び原子Bkを含む分子Bの一部において、ステップS1008で行われる原子の対の対応付けを検討する。まず、分子Aを構成する(例えば、a個の)全ての原子と、分子Bを構成する(例えば、b個の)全ての原子とから形成され得る原子同士の(a×b通りの)対のうち、原子Ai、Bkの対において、(ステップS1006で求めた)S3が、他のどの対よりも大きい、即ち最大であるとする。そうすると、まず原子Ai、Bkの対が対応付けされる。
次に、Aiを除いた分子Aを構成する(a−1)個の原子と、Bkを除いた分子Bを構成する(b−1)個の原子とから、形成され得る原子同士の(a−1)×(b−1)通りの対のうち、原子Aj、Blの対において、S3が、他のどの対よりも大きいとする。そうするとそこで原子Aj、Blの対が対応付けされる。このとき、AjはAiと直接結合しているとは限らず、BlはBjと直接結合しているとは限らない(このことは以下、同様である)。
【0047】
次に、AiとAjを除いた分子Aを構成する(a−2)個の原子と、BkとBlを除いた分子Bを構成する(b−2)個の原子とから、形成され得る原子同士の(a−2)×(b−2)通りの対のうち、原子Aj2、Bl2の対において、S3が、他のどの対よりも大きいとする。そうするとそこで原子Aj2、Bl2の対が対応付けされる。
更に次に、AiとAjとAj2を除いた分子Aを構成する(a−3)個の原子と、BkとBlとBl2を除いた分子Bを構成する(b−3)個の原子とから、形成され得る原子同士の(a−3)×(b−3)通りの対のうち、原子Aj3、Bl3の対において、S3が、他のどの対よりも大きいとする。そうするとそこで原子Aj3、Bl3の対が対応付けされる。
更に次に、AiとAjとAj2とAj3を除いた分子Aを構成する(a−4)個の原子と、BkとBlとBl2とBl3を除いた分子Bを構成する(b−4)個の原子とから、形成され得る原子同士の(a−4)×(b−4)通りの対のうち、原子Aj4、Bl4の対において、S3が、他のどの対よりも大きいとする。そうするとそこで原子Aj4、Bl4の対が対応付けされる。
【0048】
このような対応付けを、対応可能原子の対が無くなるまで、順次繰り返して行う。対応付けが終われば、その対応付けの下でのグラフマッチスコアM0を求める。
【0049】
ステップS1008では、ステップS1006で求めた多数の(例えば、a×b通りの)S3、即ちM(A、B)に基づいて、最終候補となり得る対応付け、及びその対応付けの下でのグラフマッチスコアM0の算出が行われる。
【0050】
[ステップS1010];算出したM0(A,B)が、想定され得る最大値であるか否かが確認される。具体的には、M0(A,B)が、M(A,A)又はM(B,B)に等しいかどうか、確認される。図3(2)に示すように、M(A,A)(又はM(B,B))は、最大値であると考えられるから、このステップS1010はこれ以上、グラフマッチスコアを算出する必要がないのかどうかを確認するために行われる。
【0051】
等しければ(ステップS1010・Yes)、ステップS1016にて原子対応{m(Ai)}とグラフマッチスコアM0(A、B)を出力して終了する。等しくなければ(ステップS1010・No)、ステップS1012に移行する。
【0052】
[ステップS1012];ステップS1012では、最終候補となり得る対応付けの微調整が行われる。
【0053】
分子Aにおけるひとつの結合した{Ai,Aj}の組に対応する、分子Bの{Bk,Bl}において、一方を他の原子Bnと入れ換え、入れ換えた原子についてのみ、分子Aにおける原子の対応を変更して、グラフマッチスコアM1(A,B)を求める。なお、原子Bnが分子Aにおいて対応する原子を持たない場合であってもよい。
【0054】
例えば、図9に示される、原子Ai、Aj、Amを含む分子Aの一部、及び、原子Bk、Bl、Bnを含む分子Bの一部において、ステップS1012で行われる原子の対の対応付けの変更の例を、説明する。AiとBk、AjとBl、及び、AmとBnが、対応付けられており、AiとAjが結合しているとする。ここで、{Bk、Bl}のうちの一方であるBlと、Bnとを入れ換え、AjとBnを対応付け、同時に、AmとBlを対応付ける。即ち、m(Aj)=Blであったものをm(Aj)=Bnとし、m(Am)=Bnであったものをm(Am)=Blとする。その他の原子に係る対応付けは動かされない。この一部のみ変更された対応付けに基づいて、グラフマッチスコアM1(A,B)を求める。
【0055】
図9の例における原子Bkを(図示しない)Bpと入れ替える、というような対応付けの変更であってもよい。
【0056】
[ステップS1014];算出したM1(A,B)が、M0(A,B)より大きいかどうか、確認される。即ち、ステップS1012にて、微調整を施した原子対応付けから算出されるグラフマッチスコアM1(A,B)の変動が確認される。算出したM1(A,B)が、M0(A,B)より大きければ(ステップS1014・Yes)、M1(A,B)の値がM0(A,B)に上書きされ(ステップS1015)、S1012にて更に微調整が施された原子対応付けから算出されるグラフマッチスコアM1(A,B)が求められる。
【0057】
[ステップS1016];算出したM1(A,B)が、M0(A,B)より大きくなければ(ステップS1014・No)、原子対応付けとグラフマッチスコアM0(A、B)を出力して終了する。
【0058】
5.分子の構造重ね合わせ
図4及び図2に示すフローチャートにより求めた原子対応に基づく、構造重ね合わせの表示について説明する。分子Aと分子Bの分子構造の重ね合わせにおいて、分子Aの原子{Ai}に対応した原子{m(Ai)}は適宜、重ね合わせられて表示される。このとき、Kabschの方法(McLachlan , AD. Gene duplications in the structural evolution of chymotrypsin. Journal of Molecular Biology, 128, 49-79, 1979. Kabsch, W. A solution for the best rotation to relate two sets of vectors. Acta Crystallographica, 32A, 922-923, 1976. )が用いられてもよい。
【0059】
このとき、2分子間で対応するねじれ角は、以下の方法でそろえられる。
(1) グラフマッチにより結合した分子Aの原子{Ai,Aj,Ak,Al}が、同様に結合した分子Bの原子{m(Ai),m(Aj),m(Ak),m(Al)}に対応し、かつ、結合Aj−Akと、m(Aj)−m(Ak)がいずれも回転可能結合であれば、分子Bのねじれ角{m(Ai),m(Aj),m(Ak), m(Al)}を、分子Aの対応するねじれ角{Ai,Aj,Ak,Al}と同値にする。
【0060】
上記(1)の様子を模式的に表現したのが、図6である。「回転可能結合」であるか否かの判断においては、図2におけるステップS04及びS08や、「(2.3)回転可能結合」にて設定されるデータが利用される。
【0061】
6.高速グラフマッチ探索アルゴリズムによる分子構造の重ね合わせ及びその処理の実施例。
図2に示される高速グラフマッチ探索アルゴリズムによる分子構造の重ね合わせ及びその処理のためのフローチャートを実現するプログラムを実装し、クエリ(分子A)をG39(タフミル)とし、探索ターゲットデータベースをPDBの全リガンド(9445種)として、計算を行った。動作周波数2.4GHzのデスクトップコンピュータを利用した。計算時間は、8分56秒であった。
【0062】
以下の表6で、上述の計算によりクエリの原子{Ai}に対応した{m(Ai)}を、下に並べて表示している。左端カラムには、クエリ(分子A)及び探索対象分子(分子B)の例を示している。左端から2番目のカラムには、原子数を示している。左端から3番目のカラムには、分子A(G39(タフミル))と分子Bとのグラフマッチスコアを示している。左端から4番目のカラムには、自己(自己の分子)とのグラフマッチスコアを示している。そして、右部には原子アラインメントを示している。各カラム(縦)に並んだ原子種が一致する場合は共通骨格に当たる。ここで、各カラム(縦)に並んだ原子種が90%以上一致すれば最下行に“**”を、50%以上一致すれば最下行に“++”を、示している。
【表6】

【0063】
図11は、クエリの原子クエリの原子{Ai}に対応した{m(Ai)}の組で、重ね合わせを行い原子座標を出力した図(図11(a))と、クエリ構造から共通骨格にあたる原子座標を出力した図(図11(b))である。
【0064】
7.アルゴリズム性能評価
本実施形態に係る高速グラフマッチ探索アルゴリズムの性能評価を行った。
【0065】
(7.1 アルゴリズム性能評価(1))
本実施形態は、多項式時間アルゴリズム未知のNP困難問題に近似解を与えるものである。総組み合わせ数1012以下の問題に対して、全探査により最大スコア(=正解)を求め、本実施形態に係る高速グラフマッチ探索アルゴリズムによる解と比較を行った。図12(1)は、全探査組み合わせ数に対する計算時間をグラフ化したものであり、下方から本実施形態による計算時間、全探査による計算時間、及び、本実施形態による計算時間に対する全探査による計算時間の比を示している。本実施形態に係るアルゴリズムは、グラフの探査範囲では10-4〜10-3秒で計算が可能である。全探査による場合は、10-4〜106秒を要するものである。
【0066】
図12(2)は、全探査組み合わせ数に対する正解率をグラフ化したものである。グラフの探査範囲では、平均97%の割合で正解を発見した。更に、図12(3)は、正解スコア差と累積正解率の関係をグラフ化したものである。誤答した場合でも、正解とのスコア差は最大2点であった。これら図12(1)〜(3)に示すグラフ及び数値から、本実施形態に係る高速グラフマッチ探索アルゴリズムは、高い性能を持つと考えられる。
【0067】
(7.2 アルゴリズム性能評価(2))
ブロンのクリーク探索アルゴリズム(Bron C. & Kerbosch J. Algorithm 457: Finding all cliques of an undirected graph. Communications of the Association for Computing Machinery, 16, 575-577, 1973)を用いて発見的にグラフマッチを行う方法であるsimcompの方法(Hattori, M., Okuno, Y., Goto, S. & Kanehisa, M. Development of a chemical structure comparison method for integrated analysis of chemical and genomic information in the metabolic pathways. Journal of American Chemical Society,125,11853-11865, 2003)と、成績比較を行った。
【0068】
ランダムに選んだ同じ50種の分子集合に対し総当たりグラフマッチを行い、全比較1225例(all)、及びいずれかの方法が部分グラフ(確実に正解である)を発見した136例(partial)について、本実施形態の定義によるスコアと実行時間(グラフマッチに要した実時間)を比較した。simcompの方法における最大試行回数(Rmax)を、1.5×10(デフォルト値)〜10で変化させた。
【0069】
その結果(以下、表7参照)、本実施形態に係る高速グラフマッチ探索アルゴリズムは、136の部分グラフ(partial)をすべて発見したのに対し、simcompの方法は10例(7%)で失敗した。実行時間は一例を除いて本法が高速で、平均48ミリ秒高速であった。Rmaxを増大させても発見できる部分グラフに逆転はなく、simcompの方法の実行時間が増大するだけであった。また全比較(all)においても、Rmax=1.5×10で本実施形態に係る高速グラフマッチ探索アルゴリズムが96ミリ秒遅い(但し、発見したグラフマッチのスコアは高い)以外は、どのRmaxにおいても、より高速により高スコアのグラフマッチを発見した。これらの数値から、本実施形態に係る高速グラフマッチ探索アルゴリズムは高い性能を持つと考えられる。
【表7】

【符号の説明】
【0070】
2・・・コンピュータシステム、4・・・インターネット、8・・・外部サーバ、10・・・外部データベース、12・・・出力部、14・・・中央処理部、16・・・キーボード、18・・・マウス。

【特許請求の範囲】
【請求項1】
第1の分子Aを構成する原子(Ai,Aj,・・・)の各々に係る座標データと第2の分子Bを構成する原子(Bk,Bl,・・・)の各々に係る座標データを記憶部から入力し、演算部にロードされるコンピュータプログラムに従って、演算部及び記憶部に構築される仮想メモリ空間において第1の分子Aの夫々の原子(Ai,Aj,・・・)と第2の分子Bの夫々の原子(Bk,Bl,・・・)との対応付け(m(Ai)=Bk)を求めて重ね合わせを行い(i,j,k,lはいずれも自然数)、第1の分子Aと第2の分子Bの間の最適な原子間対応、及び第1の分子Aと第2の分子Bの類似度に係るデータを出力部に出力する、第1の分子Aと第2の分子Bとの類似度を評価するための高速グラフマッチ検索装置において、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、周囲の環境が相互にどれだけ似ているかを示す第1の類似指標S1(Ai、Bk)を求める第1の算出手段と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、等しい結合距離にある周囲の原子Aj、Blの全ての組につき、第1の類似指標S1(Aj,Bl)を積算して算出する第2の類似指標S2(Ai、Bk)を求める算出手段であって、その原子Ai、Bkの対の各原子から等しい結合距離にある周囲の原子Aj、Blが同じ元素であれば、更に第1の類似指標S1(Aj,Bl)に係数を掛けた上で積算する、第2の類似指標S2(Ai、Bk)を求める第2の算出手段と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対を始点とし、第1の分子Aの原子と第2の分子Bの原子とを順次対応付けして全体の対応を作成し、そのときに算出されるグラフマッチスコアM(A,B)を値とする第3の類似指標S3(Ai、Bk)を求める算出手段であって、対応付け作成時には、既に対応付け済みの原子に直接結合する原子を次に選択すること、及び第2の類似指標S2が高い対を選択することを優先することを、条件とする、第3の類似指標S3(Ai、Bk)を求める第3の算出手段と、
第3の算出手段にて最大のS3(Ai,Bk)を算出した際の、始点の原子(Ai,Bk)の対から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の組が無くなるまで続けたときの、全体の対応におけるグラフマッチスコアM(A,B)を求める第4の算出手段と、
第4の算出手段におけるグラフマッチスコアM(A、B)が閾値より大きいならば、第1の分子Aと第2の分子Bにつき第4の算出手段で算出した原子間対応及びグラフマッチスコアM(A,B)を出力する第5の出力手段と
を含む、分子間の類似度を評価するための高速グラフマッチ検索装置。
【請求項2】
更に、
上記第4の算出手段によりグラフマッチスコアM(A,B)を算出した後に、第1の分子Aにおけるひとつの結合した{Ai,Aj}の組に対応する、第2の分子Bの{Bk,Bl}において、一方を他の原子Bnと入れ換え、入れ換えた原子についてのみ、第1の分子Aにおける原子の対応を変更して、微調整されたグラフマッチスコアM(A,B)を求める第6の算出手段を含み、
上記第6の算出手段にて算出された、微調整されたグラフマッチスコアM(A,B)が、上記第3の算出手段にて算出されたグラフマッチスコアM(A,B)より大きければ、上記第5の出力手段が、微調整されたグラフマッチスコアM(A,B)をグラフマッチスコアM(A,B)に上書きして出力を行う
ことを特徴とする請求項1に記載の高速グラフマッチ検索装置。
【請求項3】
グラフマッチスコアM(A,B)が、下記の数1で定義され、下記数1の各項は、下記の表1で定義され、下記表1の実行modeの値は入力手段を介して外部から設定されることを特徴とする
請求項1又は2に記載の高速グラフマッチ検索装置。
【数1】

【表1】

【請求項4】
上記第1の算出手段における第1の類似指標「S1(Ai,Bk)」が、下記の数2及び表2で定義され、
上記第2の算出手段における第2の類似指標「S2(Ai,Bk)」が、下記の数3及び表3で定義され、
上記第3の算出手段における第3の類似指標「S3(Ai、Bk)」が、下記の数4及び表4で定義される
ことを特徴とする請求項1乃至3のうちのいずれか一に記載の高速グラフマッチ検索装置。
【数2】

【表2】

【数3】

【表3】

【数4】

【表4】

【請求項5】
上記表3における係数が、12であることを特徴とする請求項4に記載の高速グラフマッチ検索装置。
【請求項6】
記憶部に格納される第1の分子Aを構成する原子(Ai,Aj,・・・)の各々に係る座標データと第2の分子Bを構成する原子(Bk,Bl,・・・)の各々に係る座標データを入力し、演算部にロードされる所与のコンピュータプログラムに従って、コンピュータ上に構築される仮想メモリ空間において第1の分子Aの夫々の原子(Ai,Aj,・・・)と第2の分子Bの夫々の原子(Bk,Bl,・・・)との対応付け(m(Ai)=Bk)を求めて重ね合わせを行い(i,j,k,lはいずれも自然数)、第1の分子Aと第2の分子Bの間の最適な原子間対応、及び第1の分子Aと第2の分子Bの類似度を出力部に出力する、コンピュータを用いて第1の分子Aと第2の分子Bとの類似度を評価するための高速グラフマッチ検索方法において、
第1の分子Aの全ての原子Aiと第1の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、周囲の環境が相互にどれだけ似ているかを示す第1の類似指標S1(Ai、Bk)を求める第1の工程と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、等しい結合距離にある周囲の原子Aj、Blの全ての組につき、第1の類似指標S1(Aj,Bl)を積算する第2の類似指標S2(Ai、Bk)を求める工程であって、そのAi、Bkの対の各原子から等しい結合距離にある周囲の原子Aj、Blが同じ元素であれば、更に第1の類似指標S1(Aj,Bl)に係数を掛けた上で積算する、第2の類似指標S2(Ai、Bk)を求める第2の工程と、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対を始点とし、第1の分子Aの原子と第2の分子Bの原子とを順次対応付けして全体の対応を作成し、そのときに算出されるグラフマッチスコアM(A,B)を値とする第3の類似指標S3(Ai、Bk)を求める工程であって、対応付け作成時には、既に対応付け済みの原子に直接結合する原子を次に選択すること、及び第2の類似指標S2が高い対を選択することを優先することを、条件とする、第3の類似指標S3(Ai、Bk)を求める第3の工程と、
第3の工程にて最大のS3(Ai,Bk)を算出した際の、始点の原子(Ai,Bk)の対から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の組が無くなるまで続けたときの、全体の対応におけるグラフマッチスコアM(A,B)を求める第4の工程と、
第4の工程におけるグラフマッチスコアM(A、B)が閾値より大きいならば、第1の分子Aと第2の分子Bにつき第4の工程で算出した原子間対応及びグラフマッチスコアM(A,B)を出力する第5の工程と
を含む、分子間の類似度を評価するための高速グラフマッチ検索方法。
【請求項7】
記憶部に格納される第1の分子Aを構成する原子(Ai,Aj,・・・)の各々に係る座標データと第2の分子Bを構成する原子(Bk,Bl,・・・)の各々に係る座標データを入力し、コンピュータ上に構築される仮想メモリ空間において、第1の分子Aの夫々の原子(Ai,Aj,・・・)と第2の分子Bの夫々の原子(Bk,Bl,・・・)との対応付け(m(Ai)=Bk)を求めて重ね合わせを行い(i,j,k,lはいずれも自然数)、第1の分子Aと第2の分子Bの間の最適な原子間対応、及び第1の分子Aと第2の分子Bとの類似度を評価する処理を、コンピュータに実行させるコンピュータプログラムにおいて、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、周囲の環境が相互にどれだけ似ているかを示す第1の類似指標S1(Ai、Bk)を求める第1の算出ステップと、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対の各原子からみて、等しい結合距離にある周囲の原子Aj、Blの全ての組につき、第1の類似指標S1(Aj,Bl)を積算する第2の類似指標S2(Ai、Bk)を求める算出ステップであって、そのAi、Bkの対の各原子から等しい結合距離にある周囲の原子Aj、Blが同じ元素であれば、更に第1の類似指標S1(Aj,Bl)に係数を掛けた上で積算する、第2の類似指標S2(Ai、Bk)を求める第2の算出ステップと、
第1の分子Aの全ての原子Aiと第2の分子Bの全ての原子Bkとで形成される、原子Aiと原子Bkの組の全てに関して、原子Ai、Bkの対を始点とし、第1の分子Aの原子と第2の分子Bの原子とを順次対応付けして全体の対応を作成し、そのときに算出されるグラフマッチスコアM(A,B)を値とする第3の類似指標S3(Ai、Bk)を求める算出ステップであって、対応付け作成時には、既に対応付け済みの原子に直接結合する原子を次に選択すること、及び第2の類似指標S2が高い対を選択するのを優先することを、条件とする、第3の類似指標S3(Ai、Bk)を求める第3の算出ステップと、
第3の算出ステップにて最大のS3(Ai,Bk)を算出した際の、始点の原子(Ai,Bk)の対から開始して、未対応の原子の対の中で最大のS3(Aj,Bl)を持つものを対応させることを、対応可能原子の組が無くなるまで続けたときの、全体の対応におけるグラフマッチスコアM(A,B)を求める第4の算出ステップと、
第4の算出ステップにおけるグラフマッチスコアM(A、B)が閾値より大きいならば、第1の分子Aと第2の分子Bにつき第4の工程で算出した原子間対応及びグラフマッチスコア(A,B)を出力する第5の出力ステップとを
コンピュータに実行させるコンピュータプログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate

【図8】
image rotate

【図9】
image rotate

【図10】
image rotate

【図11】
image rotate

【図12】
image rotate

【図13】
image rotate

【図14】
image rotate


【公開番号】特開2011−170444(P2011−170444A)
【公開日】平成23年9月1日(2011.9.1)
【国際特許分類】
【出願番号】特願2010−31526(P2010−31526)
【出願日】平成22年2月16日(2010.2.16)
【出願人】(503303466)学校法人関西文理総合学園 (26)
【Fターム(参考)】