NCIplot (Non-Covalent Interaction plot) の数理的基盤とアルゴリズム：密度汎関数理論からの導出と応用

4496 words

22 minutes

NCIplot (Non-Covalent Interaction plot) の数理的基盤とアルゴリズム：密度汎関数理論からの導出と応用

2026-01-04

Theoretical Chemistry

Computational Chemistry

NCIplot

Quantum Chemistry

Density Functional Theory

QTAIM

Python

last_modified: 2026-01-04

生成AIによる自動生成記事に関する免責事項: 本記事は、計算化学における非共有結合性相互作用の解析手法であるNCIplot (Non-Covalent Interaction plot) の理論的背景とアルゴリズムについて、大規模言語モデル（生成AI）によって作成された解説記事です。数式の厳密な導出や詳細な物理的解釈については、Johnson et al. (2010) および関連する原著論文を参照してください。本記事の目的は、当該手法の数学的基礎を学術的な視点から概観し、実装に向けた理解を深めることにあります。

1. 序論：電子密度に刻まれた「微弱な力」の痕跡#

化学結合の理解において、共有結合のような強く指向性のある相互作用は、分子軌道法や原子価結合理論によって明瞭に記述されてきました。一方で、水素結合、ファンデルワールス力（分散力）、立体反発、 $\pi$ - $\pi$ スタッキングといった**非共有結合性相互作用（Non-Covalent Interactions; NCI）**は、そのエネルギースケールが小さく、空間的に非局在化している場合が多いため、視覚的かつ定量的な記述が困難な課題でした。

従来、これらの相互作用の解析には、Richard F. Baderによる**分子中の原子の量子論（Quantum Theory of Atoms in Molecules; QTAIM）**が広く用いられてきました。QTAIMは電子密度のトポロジー（臨界点やボンドパス）に基づいて結合を定義する厳密な数学的枠組みですが、分散力が支配的な系や、立体障害のような「結合ではない」相互作用の可視化においては、ボンドパスが存在しない、あるいは不連続な記述になるという限界も指摘されていました。

2010年、Erin R. Johnsonらは、密度汎関数理論（DFT）で用いられる還元密度勾配（Reduced Density Gradient; RDG）という物理量を再解釈することで、これらの微弱な相互作用を実空間上で連続的な領域（Isosurface）として可視化する手法、すなわちNCI法を提唱しました。本稿では、NCI法の核となるRDGの数学的導出から立ち返り、なぜこの量が「隠れた相互作用」を炙り出すことができるのか、その数理的メカニズムを詳述します。

2. 数学的定式化：均一電子ガスからの展開#

NCI法を理解するためには、まず電子密度 $\rho(\mathbf{r})$ とその勾配 $\nabla\rho(\mathbf{r})$ が織りなす物理的景観を理解する必要があります。

2.1 Reduced Density Gradient (RDG) の導出#

密度汎関数理論（DFT）における一般化勾配近似（GGA）では、交換相関エネルギー汎関数 $E_{xc}[\rho]$ を記述するために、電子密度 $\rho$ に加えて、密度の不均一性を表す勾配 $\nabla\rho$ を導入します。ここで、無次元化された勾配として定義されるのがReduced Density Gradient (RDG, $s$ ) です。

s(\mathbf{r}) = \frac{|\nabla\rho(\mathbf{r})|}{2(3\pi^2)^{1/3}\rho(\mathbf{r})^{4/3}}

この式の分母に現れる $\rho^{4/3}$ という項は、**均一電子ガス（Uniform Electron Gas）**モデルに由来します。フェルミ波数 $k_F$ は電子密度 $\rho$ と以下の関係にあります。

k_F = (3\pi^2\rho)^{1/3}

また、均一電子ガスにおける交換エネルギー密度などのスケーリング則を考慮すると、勾配 $|\nabla\rho|$ を無次元化するための自然な尺度が $2 k_F \rho \propto \rho^{4/3}$ となることが導かれます。すなわち、 $s(\mathbf{r})$ は「局所的な電子密度の変化率が、フェルミ波長スケールと比較してどの程度大きいか」を表す指標と解釈できます。

2.2 電子密度領域における $s(\mathbf{r})$ の挙動#

RDG $s(\mathbf{r})$ の値は、分子内の異なる領域において特徴的な振る舞いを示します。

原子近傍・共有結合領域: 電子密度 $\rho$ は大きいですが、その勾配 $|\nabla\rho|$ もまた大きいため、 $s(\mathbf{r})$ は一般に有限の値を持ちますが、共有結合の中心付近（密度が極大となる点など）では $\nabla\rho \to 0$ となるため、 $s \to 0$ となることはありません（通常、共有結合上の電子密度は非常に高いため $s$ は比較的小さいがゼロにはならない、あるいは密度が高すぎてNCIの関心領域外となる）。
分子の尾部（Tail region）: 原子から十分に離れた領域では、電子密度は指数関数的に減衰します（ $\rho \propto e^{-\alpha r}$ ）。このとき、勾配は $|\nabla\rho| \propto \alpha e^{-\alpha r}$ となります。したがって、RDGは以下のように発散します。
$s \propto \frac{e^{-\alpha r}}{(e^{-\alpha r})^{4/3}} = e^{\alpha r/3} \to \infty$
つまり、電子が存在しない真空領域に近づくにつれて $s$ は増大します。
非共有結合性相互作用領域（NCI region）: ここが最も重要な領域です。ファンデルワールス力や水素結合が働く領域では、原子間の電子密度 $\rho$ は小さい値（low density）をとります。しかし、異なる原子由来の電子密度が重なり合うことで、密度勾配 $\nabla\rho$ が相殺され、ゼロに近づく（あるいは極小値をとる）領域が生じます。このとき、 $s(\mathbf{r})$ の分子 $|\nabla\rho|$ が分母 $\rho^{4/3}$ よりも速くゼロに近づくため、 $s(\mathbf{r})$ は低密度領域においてゼロに近い極小値をとります。

この「低電子密度かつ低RDG」という特異な領域こそが、非共有結合性相互作用の存在を示すシグナルとなります。Johnsonらは、縦軸に $s$ 、横軸に $\rho$ をプロット（2D NCI plot）することで、非共有結合性相互作用が $\rho \approx 0$ 近傍の鋭いスパイクとして現れることを示しました。

3. 相互作用の分類：ヘシアン行列による符号化#

RDG $s(\mathbf{r})$ がゼロに近い領域を特定するだけでは、そこに「何らかの相互作用がある」ことは分かりますが、それが引力的なのか斥力的なのかは判別できません。そこで、電子密度の2階微分である**ヘシアン行列（Hessian Matrix）**の解析が必要となります。

3.1 電子密度ヘシアン行列の固有値#

電子密度 $\rho(\mathbf{r})$ のヘシアン行列 $\mathbf{H}(\mathbf{r})$ は以下のように定義される $3 \times 3$ の対称行列です。

\mathbf{H}(\mathbf{r}) = \nabla \nabla \rho = \begin{pmatrix} \frac{\partial^2\rho}{\partial x^2} & \frac{\partial^2\rho}{\partial x \partial y} & \frac{\partial^2\rho}{\partial x \partial z} \\ \frac{\partial^2\rho}{\partial y \partial x} & \frac{\partial^2\rho}{\partial y^2} & \frac{\partial^2\rho}{\partial y \partial z} \\ \frac{\partial^2\rho}{\partial z \partial x} & \frac{\partial^2\rho}{\partial z \partial y} & \frac{\partial^2\rho}{\partial z^2} \end{pmatrix}

この行列を対角化して得られる3つの固有値を $\lambda_1 \le \lambda_2 \le \lambda_3$ とします。QTAIMの議論によれば、これらの固有値は電子密度の曲率を表しています。

3.2 第2固有値 $\lambda_2$ の物理的意味#

局所的な相互作用の性質を決定づけるのは、主に第2固有値 $\lambda_2$ です。

結合軸方向の曲率 ( $\lambda_3$ ): 原子核間を結ぶ方向では、電子密度は極小となるため、正の曲率（ $\lambda_3 > 0$ ）を持ちます。
垂直方向の曲率 ( $\lambda_1, \lambda_2$ ): 結合軸に垂直な面内での電子密度の分布を反映します。

NCI法では、 $\lambda_2$ の符号に基づいて相互作用を以下のように分類します。

引力的相互作用 (Bonding / Attractive): $\lambda_2 < 0$ 水素結合や強い分散力相互作用の場合、原子間に電子密度が蓄積（concentrate）する傾向があります。これは垂直方向の曲率が負であることを意味します。
斥力的相互作用 (Non-bonding / Repulsive): $\lambda_2 > 0$ 環状分子の中心や立体障害（Steric clash）が生じている領域では、電子密度が周囲よりも枯渇（deplete）しています。この場合、垂直方向の曲率の一つが正に転じます。
ファンデルワールス相互作用 (Van der Waals): $\lambda_2 \approx 0$ 非常に弱い相互作用では、電子密度の変化が緩やかであり、 $\lambda_2$ はゼロに近い値をとります。

3.3 $\text{sign}(\lambda_2)\rho$ による可視化#

最終的に、3次元空間内の $s(\mathbf{r})$ の等値面（Isosurface）に対し、以下の量をカラーマッピングすることで、相互作用の種類と強度を同時に可視化します。

\Omega(\mathbf{r}) = \text{sign}(\lambda_2(\mathbf{r})) \rho(\mathbf{r})

青色 ( $\Omega < 0$ , large magnitude): 強い引力（水素結合など）。
緑色 ( $\Omega \approx 0$ ): 弱い引力（ファンデルワールス力）。
赤色 ( $\Omega > 0$ , large magnitude): 強い斥力（立体障害）。

4. アルゴリズムと実装：Pythonによるプロトタイピング#

ここでは、NCIplotの計算プロセスを理解するために、Python（NumPy）を用いた簡易的な実装フローを学術的な疑似コードとして示します。実際の量子化学計算では、波動関数ファイル（.wfn, .fchk等）からグリッド上の密度情報を生成しますが、ここではグリッドデータが既に存在すると仮定した処理の流れを示します。

import numpy as np
from scipy import linalg

def calculate_nci_indices(rho_grid, grad_rho_grid, hessian_rho_grid):
    """
    NCI解析に必要な物理量を計算する関数
    
    Parameters:
    -----------
    rho_grid : ndarray (nx, ny, nz)
        3次元グリッド上の電子密度
    grad_rho_grid : ndarray (3, nx, ny, nz)
        電子密度の勾配ベクトル
    hessian_rho_grid : ndarray (3, 3, nx, ny, nz)
        電子密度のヘシアン行列

    Returns:
    --------
    rdg : ndarray
        Reduced Density Gradient (s)
    signed_rho : ndarray
        sign(lambda2) * rho
    """
    
    # 物理定数等の定義
    pi = np.pi
    # RDGの分母のスケーリング係数 2 * (3*pi^2)^(1/3)
    c_s = 2.0 * (3.0 * pi**2)**(1.0/3.0)
    
    # 1. 勾配のノルム |nabla rho| の計算
    # grad_rho_gridは (dx, dy, dz) 成分を持つ
    norm_grad = np.sqrt(np.sum(grad_rho_grid**2, axis=0))
    
    # 2. Reduced Density Gradient (s) の計算
    # rhoが極めて小さい領域でのゼロ除算を防ぐため、小さな値を加えるかマスク処理が必要
    epsilon = 1e-15
    rho_safe = np.maximum(rho_grid, epsilon)
    
    rdg = norm_grad / (c_s * rho_safe**(4.0/3.0))
    
    # 3. ヘシアン行列の固有値解析とlambda_2の抽出
    # グリッド点ごとに固有値を計算する必要がある（計算コストが高い部分）
    # ここでは概念的なベクトル化処理として記述
    
    # 配列の形状を整える (N_points, 3, 3)
    nx, ny, nz = rho_grid.shape
    hess_reshaped = hessian_rho_grid.transpose(2, 3, 4, 0, 1).reshape(-1, 3, 3)
    
    # 固有値を計算（昇順にソートされることを期待、ライブラリ依存に注意）
    # eighは対称行列・エルミート行列用
    eigvals = np.linalg.eigvalsh(hess_reshaped)
    
    # 2番目の固有値 lambda_2 を取得
    # eigvalsは (N_points, 3) で、[:, 1] が lambda_2 に相当
    lambda_2_flat = eigvals[:, 1]
    lambda_2 = lambda_2_flat.reshape(nx, ny, nz)
    
    # 4. sign(lambda_2) * rho の計算
    signed_rho = np.sign(lambda_2) * rho_grid
    
    return rdg, signed_rho

def generate_nci_isosurface(rdg, signed_rho, s_isovalue=0.5, rho_cutoff=0.05):
    """
    可視化のための等値面データのフィルタリング概念
    
    Parameters:
    -----------
    s_isovalue : float
        RDG等値面の閾値（通常 0.3 - 0.7 程度）
    rho_cutoff : float
        解析対象とする電子密度の上限（化学的に意味のある相互作用のみ抽出）
    """
    
    # 可視化ソフトウェア（VMD, ChimeraX等）では、
    # 1. RDGをスカラー場として読み込み
    # 2. RDG = s_isovalue の等値面を作成
    # 3. その等値面上に signed_rho で色付け
    # という手順を取る。
    
    # 数値解析的には、高密度領域（内殻など）や低密度すぎるノイズを除去する
    mask = (np.abs(signed_rho) < rho_cutoff)
    
    # マスク適用後のデータを返す（可視化用フォーマットへの変換は別途必要）
    return rdg * mask, signed_rho * mask

4.1 実装上の注意点#

数値安定性: 分母に $\rho^{4/3}$ があるため、電子密度が極端に小さい領域ではRDGが発散または数値的に不安定になります。適切なカットオフ（密度閾値）の設定が必須です。

グリッド解像度: NCIは微細な密度勾配の変化を捉える必要があるため、通常の電子密度プロットよりも細かいグリッド（Cube file）解像度が推奨されます。粗いグリッドでは、スパイクが消失したり、アーティファクトが生じる可能性があります。

5. 歴史的背景とインパクト：QTAIMからNCIへ#

5.1 QTAIMの功績と「隙間」#

1970年代から90年代にかけて確立されたBaderのQTAIMは、化学結合を「結合原子間の臨界点（Bond Critical Point; BCP）とそれらを結ぶボンドパス」として定義しました。これは共有結合や強いイオン結合、通常の水素結合に対しては極めて強力な記述力を持ちました。しかし、QTAIMはトポロジー（接続性）に依存するため、以下の課題がありました。

非局在化相互作用の不可視性: 例えばベンゼンのような芳香環に対する $\pi$ -スタッキングや、メタン二量体のような分散力結合では、明確なBCPが見つからないか、見つかったとしてもその物理的意味が直感的に理解しにくい場合がありました。

斥力の記述: 立体反発はエネルギー的には不安定化要因ですが、QTAIMでは「閉殻相互作用」として引力的な水素結合と同様のトポロジー分類（ $\nabla^2\rho > 0$ ）に属する場合があり、区別が困難でした。

5.2 NCIの登場とパラダイムシフト#

2010年、Yangのグループ（Johnsonら）による J. Am. Chem. Soc. 論文 “Revealing Non-Covalent Interactions” は、この状況を一変させました。彼らは、既存のDFT計算で捨てられていた情報（密度勾配）を再利用することで、計算コストをほとんど追加することなく、ファンデルワールス力や立体障害を「見て分かる」形にすることに成功しました。この手法の最大の功績は、「結合（Bond）」という二中心的な概念から、「相互作用領域（Interaction Region）」という空間的な概念への拡張です。これにより、分子認識や自己集合、タンパク質の折り畳みといった、多数の微弱な相互作用が協奏的に働く現象の解析が飛躍的に進歩しました。

6. 実利的な成果と応用例#

NCIplotの実装と普及により、以下のような化学的課題に対して具体的な知見が得られています。

6.1 立体障害の定量的可視化#

有機化学における「立体障害（Steric hindrance）」は、これまで分子模型の干渉やファンデルワールス半径の重なりとして定性的に語られてきました。NCI法は、ビシクロ系化合物や混み合った配位子を持つ遷移金属錯体において、反発領域を赤色の等値面として明示しました。これにより、反応選択性や安定性の議論に物理的な根拠を与えることが可能になりました。

6.2 水素結合ネットワークの解明#

タンパク質の二次構造（ $\alpha$ -ヘリックス、 $\beta$ -シート）や水のクラスターにおいて、水素結合は単独で存在するのではなく、ネットワークを形成します。NCI法は、これらの水素結合を連続的な青〜緑色の面として描き出し、協同効果（Cooperativity）や結合強度の分布を一目で把握することを可能にしました。また、通常の水素結合（O-H…O）と弱い水素結合（C-H…O, C-H… $\pi$ ）を同一のスケールで比較できる点も大きな利点です。6.3 実験データへの拡張近年では、高分解能X線回折実験から得られた電子密度分布に対してNCI解析を適用する試みも進んでいます（例：Saleh et al., 2012）。これにより、理論計算（気相中の孤立分子）と実験結果（結晶場中の分子）の直接比較が可能となり、結晶パッキングにおける微弱相互作用の役割解明に貢献しています。

7. 結論#

NCIplot（RDG法）は、電子密度とその勾配という基本的な物理量の中に潜んでいた「微弱な力の地図」を人類に提示しました。その数学的基盤は、均一電子ガスモデルという単純な物理モデルからの逸脱を評価するという、DFTの核心に根ざした堅牢なものです。ヘシアン行列の固有値解析を組み合わせることで、単なる密度の重なりを「引力」と「斥力」に色分けするこのアルゴリズムは、計算化学のみならず、構造生物学や材料科学においても標準的な解析ツールとしての地位を確立しています。今後、機械学習ポテンシャルとの融合や、より巨大な系への適用を通じて、その有用性はさらに拡大していくことでしょう。

参考文献#

Johnson, E. R.; Keinan, S.; Mori-Sánchez, P.; Contreras-García, J.; Cohen, A. J.; Yang, W. “Revealing Non-Covalent Interactions,” J. Am. Chem. Soc. 2010, 132, 6498–6506. (Original proposal of the RDG/NCI method)
Contreras-García, J.; Johnson, E. R.; Keinan, S.; Chaudret, R.; Piquemal, J.-P.; Beratan, D. N.; Yang, W. “NCIPLOT: A Program for Plotting Non-Covalent Interaction Regions,” J. Chem. Theory Comput. 2011, 7, 625–632. (Implementation details and the NCIPLOT program)
Bader, R. F. W. Atoms in Molecules: A Quantum Theory, Clarendon Press, Oxford, 1990. (Foundational text on QTAIM)
Saleh, G.; Gatti, C.; Lo Presti, L.; Contreras-García, J. “Revealing Non-covalent Interactions in Molecular Crystals through Their Experimental Electron Densities,” Chem. Eur. J. 2012, 18, 15523–15536. (Application to experimental electron densities)
Otero-de-la-Roza, A.; Johnson, E. R.; Contreras-García, J. “Revealing Non-Covalent Interactions in Solids: NCI Plots Revisited,” Phys. Chem. Chem. Phys. 2012, 14, 8627-8635. (Extension to solid-state systems)