今回の話題は分子の回転についてプログラミングしてみます。話が長くなりそうなので、ブログを幾つかに分割して紹介します。
はじめに
前回は母核構造を固定してコンホメーションサーチをしてみました。色々なコンホメーションを発生させた後にやってみたくなるのは母核構造の重ね合わせでしょう。例えばタンパク質結晶構造中のリガンドに、発生させた分子の母核構造を重ねてみるということなど簡単に思いつきます。思いつくのは簡単ですが、分子をうまく並進・回転させないと目的の母核構造と綺麗に重なりません。分子の並進は定ベクトルを全ての原子座標に足し算すれば良いので簡単です。しかしながら分子の回転はそんなに簡単ではありません。ましてや母核構造がピッタリ重なるように回転するのどうすればいいのでしょうか?
実を言うとこの母核構造の重ね合わせについては、OpenBabelの中にもそれらしい関数が存在することは私も確認しているのですが、APIドキュメントに一切説明がありません。使い方が良くわからないものを無理に使うよりは自分で作ろうかなと思いました。ということで、先ずは前編として分子の重ね合わせ計算のための理論を紹介します。
分子の回転
物体の3次元的な回転を扱う方法としてオイラー角と四元数の2つがあります。オイラー角の方は大抵の力学の教科書に出てくるので、物理学を学んだ事がある方は良くご存じかと思います。これは高校数学で勉強する2次元平面での回転行列$$ \begin{pmatrix} \cos \theta & -\sin \theta \\ \sin \theta & \cos \theta \end{pmatrix} $$を組みあわせて3次元の回転を表現するので、一見すると直観的な表現なのですが、回転行列には\( \sin \theta \)とか\( \cos \theta \)が複雑に出てきます。分子の重ねあわせ計算にはこれが大変面倒で、頑張れば何とかなるとは思うけど正直あまりやりたくないですね。一方、四元数の方は概念は難しいですが実際の計算の方は意外と単純です。そういうことで、ここでは四元数を使って分子の重ね合わせを考えることとします。
四元数
高校数学で勉強する複素数は実数部分を\(a\)、虚数部分を\(b\)、虚数単位を\(i\)として\(a+bi\)と表されて、虚数単位には\(i^2=-1\)という関係があります。件の四元数はこの複素数を拡張して、虚数単位として\(i, j, k\)の3種類を用います。実数部分を\(a\)、虚数部分を\(a, b, c\)とすると四元数は
$$ a+bi+cj+dk $$と表されます。そして虚数単位\(i, j, k\)の関係は以下の通りです。
$$
\begin{align}
i^2&=-1, & j^2&=-1, & k^2&=-1 \\
ij&=k, & jk&=i, & ki&=j \\
ji&=-k, & kj&=-i, & ik&=-j
\end{align}
$$
もうヤバい感じですよね(笑)。でも、面倒でも落ち着いて計算すれば大丈夫だと言うことが分かるので是非手を動かして計算してみてください。
四元数の性質
ここでは四元数の性質を簡単に紹介します。詳しくはWikipediaなどを参照してください。まず2つの四元数\(p=a+bi+cj+dk\)と\(q=w+xi+yj+zk\)の足し算です。
$$
p+q=(a+bi+cj+dk)+(w+xi+yj+zk)=(a+w)+(b+x)i+(c+y)j+(d+z)k
$$
となります。次に掛け算です。
$$
\begin{align}
pq=(a+bi+cj+dk)(w+xi+yj+zk)=aw&+axi+ayj+azk \\
&+bwi+bxi^2+byij+bzik \\
&+cwj+cxji+cyj^2+czjk \\
&+dwk+dxki+dykj+dzk^2
\end{align}
$$
となり四元数の虚数単位の関係式より以下のようになります。
$$
\begin{align}
pq=&aw-bx-cy-dz \\
&+(ax+bw+cz-dy)i \\
&+(ay-bz+cw+dx)j \\
&+(az+by-cx+dw)k
\end{align}
$$
ここで四元数の積は可換ではありません。つまり\(pq \neq qp\)です。
あと、結合法則と分配法則は成り立ちます。つまり行列と同じ代数の法則となります。
四元数をスカラー部とベクトル部に分けて表現する方法があります。虚数単位\(i,j,k\)をそれぞれ単位ベクトル\(\mathbf{i},\mathbf{j},\mathbf{k}\)とみなします。例えば上記の四元数\(p,q\)は
$$
\begin{align}
p&=a+bi+cj+dk\\
&=(a,\mathbf{p}), \qquad \mathbf{p}=(b,c,d)\\
q&=w+xi+yj+zk\\
&=(w,\mathbf{q}), \qquad \mathbf{q}=(x,y,z)\\
\end{align}
$$
と表現されます。このとき四元数の積はベクトルの内積と外積を使って
$$
pq=(aw-\mathbf{p}\cdot\mathbf{q}, a\mathbf{q}+w\mathbf{p}+\mathbf{p}\times\mathbf{q})
$$
となります。若干シンプルにまとまりましたよね。ここで、外積を使っていることから四元数の積の順番によって結果が異なることに注意してください。
共役四元数とか
複素数に複素共役があるように、四元数にも共役四元数があります。共役四元数は$$ \bar{q}=w-xi-yj-zk $$となります。四元数の掛け算\(pq\)の共役は以下のように順番が逆に入れ替わります。
$$ \overline{pq}=\bar{q}\bar{p} $$
四元数のノルムは自身の共役との掛け算となります。
$$ ||q||^2 = q\bar{q} = \bar{q}q = a^2+b^2+c^2+d^2 $$
そして四元数の掛け算のノルムは以下の通りです。
$$ ||pq||^2 = \overline{pq}pq = \bar{q}\bar{p}pq = ||p||^2\bar{q}q = ||p||^2||q||^2 $$
最後に四元数の逆元は
$$ q^{-1}=\frac{\bar{q}}{||q||^2} $$
となります。
ベクトルの回転
ノルムが1の四元数を使うとベクトルの回転を簡単に表現することができる。あるベクトル\(\mathbf{r}\)を実数部分が0の四元数として
$$ r=(0, \mathbf{r}) $$
と定義します。そしてノルムが1の四元数\(q\)をベクトル\(r\)の回転は
$$ r^{\prime}=qr\bar{q} $$
となります。ここで四元数\(q\)は三角関数を使って以下のように表現することもできて、
$$ q=(\cos\frac{\theta}{2}, \mathbf{u}\sin\frac{\theta}{2}), \quad ||\mathbf{u}||=1 $$
ここでベクトル\(\mathbf{u}\)は単位ベクトルです。このように四元数を表現した場合に上記の回転の定義式は、ベクトル\(\mathbf{u}\)を回転軸として\(\theta\)だけ右ねじの方向に回転させるという計算になります。
最後に3次元の回転行列\(R\)と四元数\(q=a+bi+cj+dk\)との関係を計算してみましょう。\(qr\bar{q}\)を頑張って計算して書き下してみると以下の式のようになります。もちろん四元数のルールである\(i^2=j^2=k^2=ijk=-1\)を適用して整理します。
$$
\begin{align}
r^{\prime}=qr\bar{q}&=(a+bi+cj+dk)(xi+yj+zk)(a-bi-cj-dk) \\
&=\left\{(a^2+b^2-c^2-d^2)x+2(bc-ad)y+2(ac+bd)z\right\}i \\
&\quad+\left\{2(ad+bc)x+(a^2-b^2+c^2-d^2)y+2(-ab+cd)z\right\}j \\
&\quad+\left\{2(bd-ac)x+2(cd+ab)y+(a^2-b^2-c^2+d^2)z\right\}k \\
\end{align}
$$
ここから行列の形にすると以下のようになります。
$$
\begin{align}
\mathbf{r}^{\prime}&=R\mathbf{r} \\
R&=
\begin{pmatrix}
a^2+b^2-c^2-d^2 & 2(bc-ad) & 2(ac+bd) \\
2(ad+bc) & a^2-b^2+c^2-d^2 & 2(-ab+cd) \\
2(bd-ac) & 2(cd+ab) & a^2-b^2-c^2+d^2
\end{pmatrix}
\end{align}
$$
とりあえず今回はここまでとします。
次回は分子を重ね合わせるための数式を導出したいと思います。
Category: OpenBabel, プログラミング関連