PIEDA Analysis of the Bromodomain in Complex with the Acetylindole Compound UZH47 (PDB ID: 6FGL)

Two things made me really happy. First, these articles(1,2) showed me how to preprocess PDB files for FMO calculations in a straightforward way. Second, thanks to faster CPUs these days, I found that FMO calculations on my MacBookAir with about 300 residues can now finish in just one day. That’s pretty amazing—20 years ago, it would have taken nearly a week using around 20 cluster machines!

An example is here

6fgl_1

6fgl_2

6fgl_3

Doing FMO Calculations on My MacBook Air: It Works!

After repeatedly rewriting input generators for FMO calculations each time I changed jobs, I decided to develop and maintain them as open-source software — and that’s how FMOkit was born.

To successfully carry out FMO calculations, we must overcome three major challenges:

  • Compiling GAMESS
  • Preparing PDB files (adding hydrogen atoms and assigning charges)
  • Generating FMO input files

In the following sections, we’ll examine each of these in detail.

Compiling GAMESS

Precompiled binaries are available for Windows and Intel-based Mac. On Linux, compilation is required, but the process is relatively straightforward. For Macs with Apple Silicon (M1–M3 chips), modifications to the source code are necessary. Please refer to the following Japanese article for detailed instructions.

Preparing PDB files (adding hydrogen atoms and assigning charges)

FMO (Fragment Molecular Orbital) calculations are a type of quantum chemical simulation, and thus require that hydrogen atoms and fragment charges be explicitly included in the molecular system. However, standard structural formats such as PDB files typically lack this information, so it must be added during preprocessing.

While commercial molecular modeling tools like Maestro or MOE are often used in pharmaceutical companies and academic labs, I opted for an open-source solution. After evaluating various options, I chose OpenMM, a molecular dynamics toolkit, as a backend for preprocessing.

You can use the script mmcifutil.py located in the utils directory to add hydrogen atoms and assign partial charges to all atoms:

python mmcifutil.py input_file output_file

For details about how the code works, please refer to this article.

Generating FMO input files

Aside from securing computational resources, generating input files for FMO calculations is arguably the most complex and unintuitive part of the entire workflow. Creating these inputs manually is not only extremely difficult, but also arguably inhumane in terms of workload and error-proneness.

Motivated by the desire to automate FMO input generation from structural data in the mmCIF format (with planned support for mol2 in the future), I began developing FMOkit.

If you prefer a graphical interface for preparing input files, I recommend using Facio, a GUI tool specifically designed for FMO workflows.

The code is shown below, using chignolin, one of the smallest known proteins, as an example.

>>> from FMOkit import System
>>> s = System(nodes=1, cores=8, memory=12000, basissets="6-31G")
>>> s.read_file("tests/5awl-addH.cif")
>>> s.prepare_fragments()
>>> with open("5awl-fmokit.inp", "w") as f:
...     f.write(s.print_fmoinput())

Although command-line interface (CLI) support is planned in the future, the code was executed from the Python console for verification purposes in this example.

The FMO calculation (RHF/MP2, 6-31G basis sets) for this 10-residue protein completed in approximately 40 minutes on my MacBook Air (Apple M3, 2024, with 16 GB of memory).

% time ~/gamess/rungms 5awl-fmokit.inp >& t5awl-fmokit.out
~/gamess/rungms 5awl-fmokit.inp >&5awl-fmokit.out  2351.21s user 142.53s system 99% cpu 41:44.37 total

gfortranでf77のfortranコードをコンパイルする

fmoutilをコンパイルする必要があったのだけど

gfortran fmoutil.f -o fmoutil

だとエラーが出るので一つオプションを追加する必要があった。

gfortran -std=legacy fmoutil.f -o fmoutil

それにしてもfortranは古のコードをよくコンパイルできるなぁと感心するのと、またコレの日々が繰り返されるのかとちょっと憂鬱ではある。

Creating mmCIF file with hydrogen atom and partial charge appended

To create an input file for FMO calculations, we need both the structural information of the protein with hydrogen atoms included and the charge information for each fragment. However, since PDB files do not contain hydrogen atoms by default, we must add them using an external tool. Additionally, because the PDB format does not provide a field for storing partial charges, we repurposed the temperature factor field for this purpose.

Since the mmCIF format allows partial charge attributes to be added without unnecessary hacks, I wrote the code in openMM and gemmi.

from openmm.app import *
from openmm import *
from openmm.unit import *
from sys import stdout
import gemmi

pdb = PDBFile("1AKI.pdb")
forcefield = ForceField('amber14-all.xml', 'amber14/tip3pfb.xml')
modeller = Modeller(pdb.topology, pdb.positions)
residues=modeller.addHydrogens(forcefield)
system = forcefield.createSystem(
    modeller.topology, 
    nonbondedMethod=PME, 
    nonbondedCutoff=1.0*nanometer, 
    constraints=HBonds)

partial_charges = []
for force in system.getForces():
    if isinstance(force, openmm.NonbondedForce):
        for atom_index in range(force.getNumParticles()):
            charge, sigma, epsilon = force.getParticleParameters(atom_index)
            partial_charges.append(str(charge)[:-2])

with open('tmp.cif', 'w') as f:
    PDBxFile.writeFile(modeller.topology, modeller.positions, f)

doc = gemmi.cif.read('tmp.cif')
block = doc.sole_block()
loop = block.find_loop('_atom_site.type_symbol').get_loop()
loop.add_columns(['_atom_site.partial_charge'], value='?')
pcs = block.find_values('_atom_site.partial_charge')
for n, x in enumerate(pcs):
    pcs[n] = partial_charges[n]
doc.write_file('output.cif')

I’m getting tired of writing FMO input generators every time I change jobs, so I’ve decided to upload them to GitHub as open-source software.

オープン x サイエンス x ケモインフォマティクス

創薬 Advent Calendar 2024の24日目の記事です。久しぶりに書いてます。

以前にアドベントカレンダーに投稿してその勢いで入門書を書きました。 おかげさまで200超えの★と90超えのforkを得られたのでまずまず貢献したのかな?と思います。

その後、CBIのケモインフォマティクスハンズオンを開催することになり、そのマテリアルもオープンにしているので興味があれば一度動かしてみてください。

ケモインフォマティクスチュートリアル 2023

一般的なヴァーチャルスクリーニング(VS)を行うための入門的なハンズオン

化合物データの前処理(脱塩等)を行い、訓練データから活性予測モデルを作成。その後VS用のライブラリを構築し、Q活性予測モデルを利用しリガンドベースのヴァーチャルスクリーニング(LBVS)を行うという内容です。

ケモインフォマティクスチュートリアル 2024

こちらは中級以上向けで、完走率は3割切るくらいだったようですが、DRYの製薬企業研究者の2,3年目が到達するレベル感だと思います(私見)。こちらは化合物の生成モデルを使っているので色々応用が効く内容になっているはずです。ワークフローはコードで書いてもKNIME/PipelienPilotにまかせてもどっちでもいいような気がしますがバージョン管理するなら前者ですかね

  1. REINVENT4の転移学習を使ってEGFR阻害剤を生成するような生成モデルを作成
  2. 生成モデルに仮想化合物群を生成させ、Gypsum-DLでSBVSの前処理をおこなう
  3. AutoDock VinaでEGFRのキナーゼドメインに対しSBVSを実行し、結果をpymolで確認する
  4. (前処理とSBVSはmaizeというflow based programmingツールを使ってワークフロー化します)

ケモインフォマティクスチュートリアル 2025

来年やるとしたらALとかHITL関連のものをやってみたいですね。この論文にあるようにエキスパートを投入して環境に手をいれるっていうAZのアイデアも面白いですし、

または、Exscientiaのように生成モデルに工夫をしてRNN+RLじゃなくてビルディングブロックと反応をうまく学習させてより現実的な構造発生させるのにAL使うっていう感じの内容もよいかなと思ってます。

それではよいクリスマスを

CBIハンズオン(chemoinformatics)のリハをおこないました

CBI年会の初日のチュートリアルセッション(TS03 ケモインフォマティクスハンズオン)のリハをおこないました。

今回は去年のように最終調整を兼ねたリハ兼ビールでも飲んで楽しむというゆるいものではなく、課題の洗い出し兼あと二週間で潰すべきタスクのリストアップに近いような作業でした。なので結構疲れたけど、このタイミングでやっておいて本当に良かったと思いました(やらんかったら、相当グダグダなセッションになっているはず)。

わたしもmac担当として、REINVENT4きちんとインストールできるように試行錯誤したり、maizeがmacでも動くようにプルリク出したりと久々に働きました。

ハンズオンの流れとしては

  1. REINVENT4の転移学習を使ってEGFR阻害剤を生成するような生成モデルを作成
  2. 生成モデルに仮想化合物群を生成させ、Gypsum-DLでSBVSの前処理をおこなう
  3. AutoDock VinaEGFRのキナーゼドメインに対しSBVSを実行し、結果をpymolで確認する
  4. (前処理とSBVSはmaizeというflow based programmingツールを使ってワークフロー化します)

という中級以上向けの攻めた内容になっていると思います。ちなみに、チュートリアルセッションは定員に達したそうで参加したい方はキャンセル待ちになるのではないでしょうか?

終了後にビールを飲みました(お約束)が、やはり来年以降はRDKitUGMJPの資金を調達しながら、ボランティアでやっているモデレーターにとってもwin-winになるような仕組みを考えたいなーと思いました。

Running REINVENT4 on the M3 chip

As it could not be installed with the original requirements, I modified it and uploaded it to the Mishima-syk repository.

$ git clone https://github.com/Mishima-syk/REINVENT4.git
$ cd REINVENT4
$ conda create --name reinvent4 python=3.11
$ conda activate reinvent4
$ pip install -r requirements-macOS.lock
$ pip install --no-deps .

Now, you can set mps (Metal Performance Shaders) in the device parameter

# mps_sampling.toml
# REINVENT4 TOML input example for sampling
run_type = "sampling"
device = "mps"  # M3 GPU
json_out_config = "_sampling.json"  # write this TOML to JSON

[parameters]

## Reinvent: de novo sampling
model_file = "priors/reinvent.prior"
output_file = 'sampling.csv'  # sampled SMILES and NLL in CSV format
num_smiles = 157  # number of SMILES to be sampled, 1 per input SMILES
unique_molecules = true  # if true remove all duplicatesd canonicalize smiles
randomize_smiles = true # if true shuffle atoms in SMILES randomly

Running the reinvent4 program

$ reinvent -l sampling.log mps_sampling.toml

Mishima.syk #21やります

ちょっと間が空いてしまいましたが、MIshima.syk #21をやります。

発表したい方は早めに時間を取っておいてください(いつも駆け込みでぶっこんでくるので、、、) それからオンラインで参加したい場合には私に連絡をください。zoomのアドレスを送ります。

今回はLLMの話とか、ケモインフォマティクスの深い話とかになりそうですが、今回もディープな話で盛り上がれることを期待しています。

余談ですが、今JCUP XIIに出ているのですが、知り合いからラボの先輩がYouTuberになって量子化学に関して熱く語っていという話を聞いて驚きました。

そして私の講演は難しくてわからんとも。わかりやすいように周辺知識も追加したのに、、、 ちなみにトークでの私の主張は以下です。

  • Our aim is to describe the ligand association/dissociation process as a chemical reaction, where reactants and products remain unchanged.
  • LBDD can be regarded as SBDD in ligand-protein complex is unknown.

共感してくれる人がいると嬉しいのですが。

MacにChEMBLを入れる

ちょっとChEMBLを使う必要ができたので、久しぶりに手元のMacbookにインストールしました。最新のバージョンは34です。

Postgresqlのインストール

バージョンが14と古かったので最新の16を入れた。

brew install postgresql@16
echo 'export PATH="/usr/local/opt/postgresql@16/bin:$PATH"' >> ~/.zshrc #psqlがない場合
source ~/.zshrc
brew services start postgresql@16
brew services list #動いているか確認

動いていたらChEMBLのデータを入れる

createdb chembl_34
pg_restore --no-owner -d chembl_34 chembl_34_postgresql.dmp #30分くらいかかる

これでOK

DBにアクセスするには

psql chembl_34

「機械学習による分子最適化」を読んだ

京都に行く新幹線の中で機械学習による分子最適化を読みました。

しっかりと数式で説明されていて簡単に概要を理解するというよりは、ある程度このあたりわかっている人が現場を制するための書籍のように思いました。大体PRMLは必須で、ECFPのアルゴリズムも知っていることが前提、分子生成系の論文をある程度読んでいる方が好ましいといった感じ。あとは強化学習も読んでおきましょう。

SELFIESを丁寧に解説していて個人的には参考になりました。CFGとかlisp処理系作って遊んでたときに勉強したけど、ここでも使われてんのかーとか。ただ、実践的にはSMILESで十分な気もする。

ニューラルグラフフィンガープリントは改めてハッシュ関数をニューラルネットワークで置き換えることで学習可能なパラメタを与えつつ微分可能っていう記述を目にして、ちょっと面白いパラメタの与え方を思いついたので今度ためしてみる。

それからオンライン強化学習とオフライン強化学習に関しても改めて自分の中で整理された。製薬企業のケモインフォマティクス研究のガチの人が読むと面白いかと思います。というか、このあたりサラサラ読まんとガチとは言えんよね的な本ですかねぇ。

ただ、自分の中ではこっちで頑張るって決めてるんで、なんか精度が明らかに上がるようなコンセプトを提案できると嬉しいなぁと思っています。