速報 (Selected Paper)
特徴ベクトルから分子構造を生成するアルゴリズム
キーワード:
KEYWORDS Machine learning,
Feature vector,
Molecular graph reconstruction,
De Bruijn graph,
Lattice basis
2018 年 17 巻 5 号 p. 228-231
詳細
2018 年 17 巻 5 号 p. 228-231
Machine learning for predicting properties of molecules and designing molecules with desired properties requires suitable numerical representation of molecular structure. In this study we propose to represent a molecular graph as a high dimensional feature vector by counting the number of substructures of a molecule, and report conditions that these vectors should satisfy. The vectors are expanded by the lattice basis extracted from the feature vectors of a large molecular dataset. In addition, an efficient algorithm to reconstruct molecular structures from the feature vector is presented. The latent space provides a continuous and concise representation of molecular graphs, and the representation allows smooth interpolation between two molecules. The proposed method has wide application to studies of quantitative structure activity relationships and machine learning of molecules.
分子構造から活性を予測したり，予測モデルから高活性の分子構造を推定する方法として，構造活性相関が古くから知られている．最近，機械学習で画像などを学習し分類したり，学習結果から新しい画像を生成する研究が実用化しつつある．また機械学習を構造活性相関に用いて，構造から活性を導くモデルが提案されている．しかし活性から構造を生成するには難しい課題がある．
SMILES文字列で分子を表し，それをneural network (NN) で学習して，分子構造を潜在変数の空間で表す報告がある [1]．この方法の難点は，似た分子がかなり異なるSMILESを与える場合，分子構造と潜在変数の間の写像が不連続になる事，また潜在変数からSMILESを生成すると，分子に対応しないものがしばしば生成される事である．分子構造を頂点(原子)と辺(結合)で構成されるグラフで表現し，頂点や辺の出現確率を学習すると，生成モデルが得られる [2]．この方法では，分子を辺の重みつき完全結合グラフとして表すため，大きな分子が生成できず，また無効分子も生成される事が問題である．
Faulonらは以前，拡張原子価に基づいたsignatureと呼ばれる新しい記述子や構造活性相関モデルを提示した [3]．Signatureは予め選んだ数hに対して，各原子のh近傍原子までの部分構造の数を並べたベクトルで分子を表現する．これは分子グラフ構造から性質を予測するWeisfeiler-Lehman graph kernel [4] やmessage passing NN [5] と等価な入力に相当する．最小のh=1の場合でもsignatureは多種類あるので，この特徴ベクトルは高次元となる．このため表現が分子構造に対応する条件や，表現から分子構造を再構成する事が困難だった．我々は，表現が分子構造を持つ必要十分条件と，分子構造を再構成する新しい効率的アルゴリズムを提案する．
我々は文献 [3]を参考に特徴ベクトルの部分木構造を設計した．1,3-cyclopentadione分子のh=1の部分木構造の例をFigure 1に示す．Figure 1 (a)が示すように，分子構造の水素原子を省略し，各原子には結合した非水素原子数を追加する．部分木(b)-(e)の全原子にもこの数が書かれ，また結合次数も記される点が，先行研究と異なる．この違いが，特徴ベクトルが分子構造を持つ条件や，特徴ベクトルから分子構造を再構成する事を簡単にする．
Examples of h=1 subtrees of molecule. (a) shows the molecular graph of 1,3-cyclopentadione. A node of it is composed of an element symbol and the number of bonded non-hydrogen atoms. An edge describes the bond type. All the subtrees except for (e) appeared twice in the molecule.
一般にsignature などの特徴ベクトルは，対応する原像(分子構造)を持つとは限らない．Signature が原像を持つ必要条件は報告されている [6]．Figure 1 の分子ではその条件は，C3=O1結合に注目すると(b)と(c)の構造は同数あり，C2-C2結合に注目すると(d)の構造が偶数個になる事などである．逆に全種類の結合がこの条件を満たす時，対となる結合を融合して，少なくとも1個のグラフが得られる．しかしこのグラフは分子構造としては存在しない自己ループや多重辺を含む．我々は数種類の部分グラフを禁止する事で，特徴ベクトルが原像を持つ必要十分条件を得た．
我々のh=1の特徴ベクトルはC, N, O, ハロゲンなどの典型元素を対象とすると，1000 次元程度の非負整数ベクトルとなる．原像を持つ必要条件は，要素(部分構造数)間に多数の線形等式条件を課す．つまり原像がある解を求めるには，この条件を表す行列
右辺ベクトル
分子A, Bの特徴ベクトル
次に，特徴ベクトルが指示する部分構造を組み合わせ，分子グラフを生成するアルゴリズムを考案した．ショットガン法でDNA配列を決定する際，数塩基からなる断片の重なりをde Bruijn グラフにまとめ，それの一筆書き経路(Euler閉路)から，全体配列を決定する [8]．我々は部分構造の重なりを同様のグラフにまとめ，全てのEuler閉路を求め，分子構造を作るプログラムを作成した．そのアルゴリズムは次のようになる．
分子構造再構成アルゴリズム：
(1) 各部分木構造のroot atomが奇の次数なら，dummy atom Xiを付け，Xiをroot atomとする木を追加する．
(2) 4次以上の偶の次数のroot atomを複製し分割し，次数を全て2とする．複製相手を記録する．これらの部分構造は分子中の3原子の結合を表す．
(3) これら部分構造のde Bruijn グラフを作る．
(4) de Bruijn グラフの全てのEuler閉路を求める．
(5) 求めた閉路の辺はroot atomを示す．隣接root atomを結合しグラフを作る．
(6) 複製された次数が4以上のroot atomを融合し，dummy atomを消去する．
Figure 2の例で説明する．Figure 2 (d)の分子は，Figure 2 (a)黒で示した部分構造を持つ．この内第2, 3番目の構造は2個ある．以下特徴ベクトルは原像を持つ必要条件(*)を満たすとする．まずステップ(1)でグラフの全頂点を偶次数とし，Euler閉路を持たせる．仮定(*)からこれら部分構造をまとめたグラフは存在する．そこで握手定理より奇次数のroot atomは偶数個存在し，ステップ(1)は必ず成功する．この操作で赤で記した構造が加わり，Figure 2 (d)はEuler閉路を持つ．違うdummy atomの入れ方では，違うEuler閉路を持ち得るが，最終的に同じ分子が得られる事が分かる．ステップ(2)では次数4のroot atomを持つ5番目の部分構造を，次数2の2個の部分構造に分ける．この分け方(3通り)は，Figure 2 (d)のEuler閉路でN3の通り方を指定するが，どれも最終的に同じ分子が得られる事が分かる．ステップ(3)では，まず第1番目の部分構造と連結できる任意の構造を選び，de Brujin グラフを作る．仮定(*)よりこのような構造は常に見つかり，全部分構造をまとめたFigure 2 (b) のde Brujin グラフができる．Figure 2 (c)はそのEuler 閉路を示す．ステップ(5)と(6)によりFigure 2 (d)の分子を表すグラフが得られる．
Description of the reconstruction algorithm. Figure (a) shows a set of subtrees after step (1). From these subtrees steps (2) and (3) yield the de Bruijn graph as shown in (b). Figure (c) is a Euler cycle of it, and through steps (5) and (6) we obtain the reconstructed molecular structure (d) with a dummy node X1.
我々はtyrosine-protein kinase (JAK2)を標的とする 分子をCHEMBL [9]で検索し，その内ALOGPとIC50が記録され，かつ矛盾した実験結果がない4013分子を選んだ．またGDB13 データセット [10]から最大7個の非水素原子を持つ7211分子を選び，11224分子を得た．これら分子のh=1の特徴ベクトルを作り，LLLアルゴリズム [7, 11]で格子基底を計算した．特徴ベクトルは1099次元で，簡約により938本の格子基底を得た．つまり11224分子の特徴ベクトルは，938本の格子基底の整数結合で表現でき，11224分子を938次元の一点に写像した事になる．この計算にはIntel Core i7-5930K CPU (6 cores, 3.30GHz)を用い6時間程度で終了した．得られた基底は理解難いものもあるが，意味の分かり易いものも存在した．例えば20番の基底
次に二分子に対応する格子点間を渡り歩き，分子の内挿を試みた．提案アルゴリズムを用いて内挿した特徴ベクトルから分子構造を再構成する時，複数構造が得られる場合，似た骨格構造のものを選んだ．その結果をFigure 3に示す．特徴ベクトルに10本の格子基底を順に加えると，左下分子が右上分子に徐々に変化する様子が分かる．
Interpolation of molecules in the feature vector space. The feature vector of the lower left molecule is transformed to that of the upper right one by adding ten lattice bases,
以上に示したように，我々は適切に設計した特徴ベクトルを用いて，分子を部分木構造の集まりとして表し，高次元格子点に写像した．特徴ベクトルが原像を持つ条件を求め，また格子点から分子構造を再構成する効率的なアルゴリズムを提示した．分子の活性を特徴ベクトルと回帰する事で，最適な分子設計は整数計画問題に還元され，提案アルゴリズムで効率的に解けると期待される．
