magattacaのブログ

日付以外誤報

RDKitで3次元構造に挑戦! part 3  〜HALAVENの共結晶構造との比較〜

前回までにRDKitを用いた3次元構造の生成、ならびに2次元構造の類似性比較を行いました。
今回はPDBに登録された共結晶構造との3次元構造の比較を行いたいと思います。

1. HALAVEN のおさらい

まずは今更ですが、そもそもHALAVEN®️( eribulin mesylate)ってどういう医薬品なの??という話から。
以下インタビューフォーム(IF)より抜粋
①ハラヴェンの有効成分エリブリンメシル酸塩はHalichondrin Bの合成誘導体 。(IF p.1)
②再発乳癌等を適応とする抗悪性腫瘍剤。(IF p.11)
③作用機序としては微小管のプラス端(伸長端)に結合することによる、微小管ダイナミクスの阻害。
 微小管による正常な紡錘体形成を妨げることで、細胞分裂を停止させアポトーシスによる細胞死を誘導し、がん細胞の増殖を抑制します。(IF p.23)
④薬効評価の in vitro 試験(無細胞系)として、微小管ダイナミクスに対する作用が評価されており、微小管を形成するタンパク質であるチューブリンの重合阻害試験が用いられています。(IF p.25)

ちょうど良いことにチューブリン(Tublin) とHALAVEN®︎(Eribulin)の共結晶構造[PDB ID: 5JH7]が報告されていますので、早速眺めてみましょう。

(チューブリン命名したの日本人研究者の方だそうです。「チューブラリン」と迷ったってホンマかいな・・・チューブリン - Wikipedia)

2. 共結晶構造を眺める

まずはどのような構造となっているかPDB viewerで確認してみましょう。

f:id:magattaca:20181125151021p:plain
Fig. 1 Tublin-Eribulin 共結晶構造 (PDB ID : 5JH7)

図のようにαβ-tublinの複合体の間にEribulinが挟まった構造となっています。
Eribulinはちょうど「お椀」のような形となっており、β-tublin をキャップするかのように結合し、重合の伸長を阻害している様がわかります。
念のためcootでEribulin結合部位の電子密度を確認して見ましたが、リガンドと重なるよう電子密度を確認できました。

また、こちらの共結晶構造にはGTP/GDPも含まれており、チューブリンの伸長がGTPの結合・加水分解により駆動されている様もわかる興味深い構造となっています。

3. RDKitでリガンド構造をながめる

それではPDBに登録されているリガンドの構造をRDKitで眺めて見ましょう。
RSCB PDB のページより[Ligand ID : 6K9]のsdfファイルをダウンロードしました。

「Ligands_noHydrogens_withMissing_2_Instances.sdf」 という不穏な名前で保存されました・・・
モヤモヤしますがとりあえず「6K9_noH.sdf」とファイル名を変更して進めていきます。

まずはSDFファイルの読み込みから・・・(参考記事①)

#読み込み方法-1  (SDMolSupplier)
from rdkit import Chem
suppl = Chem.SDMolSupplier('./6K9_noH.sdf')
mols = [x for x in suppl if x is not None]  #分子が読み込めているかを確認
len(mols)  #2と出力。2つ入っているらしい

#読み込み方法-2  (Pandasデータフレームへの読み込み)
from rdkit.Chem import PandasTools 
df = PandasTools.LoadSDF('./6K9_noH.sdf')
#MolオブジェクトはROMolというカラムに格納されているらしい
eri_PDB_noH = df.loc[0, 'ROMol']

#2次元描画
from rdkit.Chem import rdDepictor
from rdkit.Chem.Draw import IPythonConsole
rdDepictor.SetPreferCoordGen(True)
eri_PDB_noH

#3次元描画
IPythonConsole.drawMol3D(eri_PDB_noH)

以下のような構造が出力されました。

f:id:magattaca:20181125162318p:plain
Fig.2 PDBに登録されているEribulin (ID: 6K9) の構造

以前、PDBフォーマットの見方で確認したように、PDBの構造には水素原子が含まれていません。
Pandasデータフレームをみると、リガンドのsdfファイルには不斉や結合次数の情報は含まれているようですがやはり水素原子は含まれていないようです。
水素原子を付加してみましょう。

#水素原子を付加
from rdkit.Chem import AllChem
eri_PDB_addH = AllChem.AddHs(eri_PDB_noH)
IPythonConsole.drawMol3D(eri_PDB_addH)

f:id:magattaca:20181125163846p:plain
Fig.3 水素原子の付加を試みた

・・・疾走感 !!!
どうやら3次元構造全体に対して水素原子1つを付け足しただけのようです。
AddHsの使い方を間違えたらしい。

いろいろ探しているといままで気づいていなかったページがPDBにありました。

f:id:magattaca:20181127004724p:plain
Fig.4 RSCB PDB リガンド情報について

こちらのリガンド詳細情報のページには
・Structure Data File (Ideal SDF)
・Structure Data File (Model SDF)
という、2つのSDF形式の情報があり、共に水素原子を含むものでした。
それぞれを
・6K9_SDF_ideal.sdf
・6K9_SDF_model.sdf
という名前で保存して、立体構造をながめてみましょう。
また、折角なのでRDKitでsmilesから生成させた3次元構造もならべて見ましょう。

#水素原子が大事なので 「removeHs=False」として読み込む
SDF_ideal = Chem.SDMolSupplier('./6K9_SDF_ideal.sdf', removeHs=False)
SDF_model = Chem.SDMolSupplier('./6K9_SDF_model.sdf', removeHs=False)

#それぞれMolBlock形式に変換
SDF_ideal_mb = Chem.MolToMolBlock(SDF_ideal[0])
SDF_model_mb = Chem.MolToMolBlock(SDF_model[0])

#smilesからの3次元構造生成は省略 ( eri_molH_mb という名前でMolBlockを生成した)
#py3Dmolで3つ並べて表示
import py3Dmol

v = py3Dmol.view(width=1000, height=400, linked=True, viewergrid=(1,3))
v.addModel(eri_molH_mb, viewer=(0, 0))
v.addModel(SDF_ideal_mb, viewer=(0, 1))
v.addModel(SDF_model_mb, viewer=(0, 2))
v.setBackgroundColor('0xeeeeee')
v.setStyle({'stick':{}})
v.zoomTo()
v.show()

以下のようになりました。
左から①RDKitで生成した構造、②PDBのIdeal SDF、③ PDBのModel SDF となります。

f:id:magattaca:20181125174853p:plain
Fig. 5 様々なソース由来の3次元構造を比較


Fig. 5 より①と②が類似しており、③は水素原子のないPDB構造(Fig. 2)に似ていることがわかります。
おそらく、Fig. 5 ②は2次元構造式から計算により生成された3次元構造で、
Fig. 5 ③ が Fig. 2の実験由来のモデル構造に水素原子を付加した形となっているのではないでしょうか?(すみません。憶測です。)

折角なので類似の①と②を重ね合わせて確認して見ましょう。

#水素原子を含まないPDBのSDFファイルをMolBlock形式に変換
eri_PDB_noH_mb = Chem.MolToMolBlock(eri_PDB_noH)

#Model SDFの構造に上記を重ね合わせる
AllChem.AlignMol(eri_PDB_noH, SDF_model[0])

#一つのviewオブジェクトに表示
p = py3Dmol.view(width=400,height=400)
p.addModel(eri_PDB_noH_mb)
p.addModel(SDF_model_mb)

p.setStyle({'stick':{}})
p.setBackgroundColor('0xeeeeee')
p.zoomTo()
p.show()

f:id:magattaca:20181125181226p:plain
Fig. 6 Model SDFとPDB ligandのSDFファイルの重ね合わせ

・・・完全に一致!
Model SDFが実験データに水素原子の情報を付け足したものと考えて良さそうです。

ところでRDKit で生成した3次元構造(Fig.5 ①) は PDBのリガンド構造(Fig. 2)よりも平面的となっており、
計算結果と実験結果との差異が見られて興味深い結果となりました。
Fig. 1 から、EribulinはTublin をキャップするように重合阻害している様子でしたので、計算で出した構造ではあまりそのような構造の特徴が再現できていないように見えます。

4. 立体配座生成のアルゴリズムを比較

念のため、他のアルゴリズムでも3次元構造を生成して比較して見ましょう。(参考記事②)

#isomeric SMILES取得過程は省略(eri_isoという変数に格納)
#それぞれ3次元構造を生成しMolBlock形式に変換

#①ディスタンスジオメトリー(DG)法 + Universal Force Field(UFF)  
eribulin_DG = AllChem.AddHs(Chem.MolFromSmiles(eri_iso))
AllChem.EmbedMolecule(eribulin_DG)
AllChem.UFFHasAllMoleculeParams(eribulin_DG) #True
AllChem.UFFOptimizeMolecule(eribulin_DG_UFF)
eribulin_DG_UFF_mb = Chem.MolToMolBlock(eribulin_DG_UFF)

#②ディスタンスジオメトリー(DG)法 + Merck Molecular Force Field (MMFF) 94 
eribulin_DG_MMFF = AllChem.AddHs(Chem.MolFromSmiles(eri_iso))
AllChem.EmbedMolecule(eribulin_DG_MMFF)
AllChem.MMFFOptimizeMolecule(eribulin_DG_MMFF)
eribulin_DG_MMFF_mb = Chem.MolToMolBlock(eribulin_DG_MMFF)

#③KDG (DG + basic Knowledge)
eribulin_KDG = AllChem.AddHs(Chem.MolFromSmiles(eri_iso))
AllChem.EmbedMolecule(eribulin_KDG, AllChem.KDG())
eribulin_KDG_mb = Chem.MolToMolBlock(eribulin_KDG)

#④ETDG (DG + Experimental Torsion angle preferences)
eribulin_ETDG = AllChem.AddHs(Chem.MolFromSmiles(eri_iso))
AllChem.EmbedMolecule(eribulin_ETDG, AllChem.ETDG())
eribulin_ETDG_mb = Chem.MolToMolBlock(eribulin_ETDG)

#⑤ETKDG (DG + K + ET)
eribulin_ETKDG = AllChem.AddHs(Chem.MolFromSmiles(eri_iso))
AllChem.EmbedMolecule(eribulin_ETKDG, AllChem.ETKDG())
eribulin_ETKDG_mb = Chem.MolToMolBlock(eribulin_ETKDG)

#py3Dmolで描画(省略)

#PDBから取得した構造(Model SDF)とのRMSDを計算
mol_list = [eribulin_DG_UFF, eribulin_DG_MMFF, eribulin_KDG, eribulin_ETDG, eribulin_ETKDG]
RMSD = []
for i in mol_list:
    a = AllChem.GetBestRMS(SDF_model[0], i)
    RMSD.append(a)

#水素原子がない場合(Heavt Atom RMSD)も計算
Heavy_Atom_RMSD = []
for i in mol_list:
    a = AllChem.GetBestRMS(eri_PDB_noH, Chem.RemoveHs(i))
    Heavy_Atom_RMSD.append(a)

f:id:magattaca:20181126222708p:plain
Fig.7 異なるアルゴリズムで生成した構造とPDB共結晶構造との比較


今回の結果では⑤ETKDG が最もRMSDが小さくなりました。
DG +MMFF(②) > KDG(③) > DG + UFF(①) > ETDG(④) > ETKDG(⑤) となっています。
概ね改良された順番にPDBの共結晶構造に近い3次元構造を生成しているようです。

最もRMSDが大きい②と、最も小さい⑤についてそれぞれPDBの構造をと重ね合わせて見ましょう。

#Model SDFの構造に②DG+MMFF と ⑤ETKDGをそれぞれ重ね合わせる
AllChem.AlignMol(eribulin_DG_MMFF, SDF_model[0])
AllChem.AlignMol(eribulin_ETKDG, SDF_model[0])

#MolBlock形式に変換
eribulin_DG_MMFF_mb = Chem.MolToMolBlock(eribulin_DG_MMFF)
eribulin_ETKDG_mb = Chem.MolToMolBlock(eribulin_ETKDG)
SDF_model_mb = Chem.MolToMolBlock(SDF_model[0])


#描画 
v = py3Dmol.view(width=1000, height=400, linked=True, viewergrid=(1,2))
v.addModel(eribulin_DG_MMFF_mb, viewer=(0, 0))
v.addModel(SDF_model_mb, viewer=(0, 0))

v.addModel(eribulin_ETKDG_mb, viewer=(0, 1))
v.addModel(SDF_model_mb, viewer=(0, 1))

v.setBackgroundColor('0xeeeeee')
v.setStyle({'stick':{}})
v.zoomTo()
v.show()

こんな感じ

f:id:magattaca:20181126224317p:plain
Fig. 8 RMSD 最大と最小をそれぞれ重ねあわせ

左右を比較してみるとRMSDの小さい右の方が重なりが良いというのがはっきりとわかります。

「化学の新しいカタチ」さんの記事(参考記事②)でDr Tomさんが 述べられているように、まずは⑤ETKDGを試しておけば良い結果が得られそうです。

5.EribulinとPaclitaxelの比較

これまでEribulin や元になった天然物 Halichondrin Bとの比較を行ってきました。
最後に「化学の新しいカタチ」さんでも何度か取り上げられている、複雑な天然物の代名詞 Paclitaxel とEribulin を比較して見たいと思います。

興味深いことにPaclitaxel(Taxol)も微小管を構成するチューブリンのβ-サブユニットに結合し、細胞分裂を阻害することで抗腫瘍活性を発揮します。
一方で、Eribulinは微小管の重合を阻害するのに対し、Paclitaxelは微小管を安定化させ脱重合を阻害するという、メカニズム上の違いもあります。

PDBを検索しましたが、残念ながらX線結晶構造は見当たりませんでしたので、Cryo-EMのうち比較的新しい[PDB id : 6EW0]を眺めて見ます。

f:id:magattaca:20181126233946p:plain
Fig. 9 微小管に結合するPaclitaxel(Taxol) [6EW0]

Eliburinはβ-Tublinに対して伸長軸をキャップして妨げるように結合していましたが、Paclitaxelは伸長軸に沿って側面から結合しています。

2次元および3次元構造を比較しましょう。
3次元構造はいずれもPDBから水素原子の付加されたmodelのsdfファイルを取得して描画しています。

f:id:magattaca:20181126235417p:plain
Fig. 10 Eribulin(左)、Paclitaxel(右)

いずれも全合成が成し遂げられているのが不思議でたまらないほど複雑な構造をとっています。
この両者の類似性をFingerprintで評価したらどのような結果になるのでしょうか??

Morgan Fingerprint (radius =2, ビット長 =2048) で類似性を評価してみます。

f:id:magattaca:20181127001822p:plain
Fig. 11 Morgan Fingerprint による原子ごとの寄与率を可視化

案の定、Tanimoto係数、Dice係数ともに 0.1 ~ 0.2 と非常に低い値となっています。

・・・まあ標的は同じタンパク質でも結合している場所全然違いますし

原子ごとの寄与率の描画では複雑な縮環部位の寄与が大きくみえるのが面白いですね。

まとめ

以上、3回にわたってRDKitをもちいた構造の描画、立体構造の取り扱いを勉強してみました。
PDB viewerでは取り扱えない、低分子リガンド構造を様々な視点から検証でき、まだまだいろいろ遊べそうです。

・・・結局コピペしかしてないので自分ではコードの書き方がわからない

jupyter notebookに打ち込みつつ記載しましたが、使い方や解釈の仕方でいろいろと間違っている部分があるかと思います。ご指摘いただければ幸いです。


先ほど向山先生の訃報を目にしました。ご冥福をお祈りします。
直接の面識はもちろんありませんが、私のような有機化学の末端員でも数々の伝説を耳にしました。
向山先生によるTaxolの全合成は英語版のWikipediaに詳しいです。
Mukaiyama Taxol total synthesis - Wikipedia

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

(1) ハラヴェンに関する情報
医薬品インタビューフォーム (pmdaのページ)
ハラヴェン静注1mg
KEGG
医療用医薬品 : ハラヴェン


(2) 共結晶構造
Tubulin-Eribulin complex [PDB ID: 5JH7]
Eribulin [Ligand ID : 6K9]
RCSB PDB - 5JH7: Tubulin-Eribulin complex
元の文献(open access
Termination of Protfilament Elongation by Eribulin INduces Lattice Defects that Promote Microtubule Catastrophes.
Doodhi, H., et. al., (2016) Curr.Biol. 26: 1713-1721
https://www.cell.com/current-biology/fulltext/S0960-9822(16)30407-9?_returnURL=https%3A%2F%2Flinkinghub.elsevier.com%2Fretrieve%2Fpii%2FS0960982216304079%3Fshowall%3Dtrue


(3)「化学の新しいカタチ」さんの以下の記事を参考にさせていただきました。
参考記事① : RDKitでケモインフォマティクスに入門
参考記事② : RDKitによる3次元構造の生成


(4) RDKit ブログ / Documentation
RDKit: 3D structure of a natural product
Getting Started with the RDKit in Python — The RDKit 2017.09.1.dev1 documentation

(5) RCSB PDB
PDB ID: 5JH7
Termination of Protofilament Elongation by Eribulin Induces Lattice Defects that Promote Microtubule Catastrophes.
Doodhi, H., et. al., (2016) Curr.Biol. 26: 1713-1721

PDB ID:6EW0
Cryo-EM structure of GDP-microtubule co-polymerised with doublecortin and supplemented with Taxol
Manka, S.W., Moores, C.A. (2018) Nat. Struct. Mol. Biol. 25: 607-615

NGL Viewer
AS Rose et al. (2018) NGL viewer: web-based molecular graphics for large complexes.
Bioinformatics dio:10.1093/bioinformatics/bty419


~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~