Zenn
⚗️

逆合成ちゃんを支える技術

2020/12/14に公開
Python
Twitter
逆合成解析
cheminformatics
tech

これは 創薬 (dry) Advent Calendar 2020 の15日目の記事です。

この記事では先日公開した逆合成ちゃんを支える技術について解説していきます。

ちなみに、逆合成ちゃんは #souyakuchan_magic から着想を得て生まれました。創薬ちゃんさんにこの場を借りてお礼申し上げます。

AiZynthFinderの使い方

逆合成ちゃんはオープンソースの逆合成ツールAiZynthFinderを手軽に試せるようにするために作りました。AiZynthFinderの詳細については別の記事で後日解説する予定なので、ここでは手っ取り早く動かす方法について紹介します。

基本的にREADMEを読めば済む話なのですが、いくつかハマったポイントがあるので補足しつつ解説していきます。

インストール

AiZynthFinderはRDKit, TensorFlow, Graphvizに依存しています。RDKitはpipでインストールできない（RDKitのブログ Why the RDKit isn't available on PyPi を参照）ので、まず何らかの方法でAnacondaをインストールする必要があります。

Anacondaはどうインストールしてもいいですが、ここではMiniconda 3の最新版をインストールすることにします。筆者はLinuxを使っているので

wget https://repo.anaconda.com/miniconda/Miniconda3-latest-Linux-x86_64.sh
bash Miniconda3-latest-Linux-x86_64.sh

のようにしました。

ハマりポイント1: Pythonのバージョン

AiZynthFinderのREADMEにはPython 3.6, 3.7 or 3.8に対応していると書いてありますが、2020年12月現在、RDKitの最新版はPython 3.8に対応していません。Anacondaの最新版ではPython 3.8がインストールされるのでPythonのバージョンを下げる必要がありました。Anacondaで使うPythonのバージョンを3.7にするには次のコマンドを実行します。

conda install python=3.7

依存ライブラリは次のコマンドでインストールできます。GPUがついているマシンではtensorflow-gpuを指定してください。

conda install -c rdkit "rdkit=>2019.09.1" -y
conda install -c anaconda "tensorflow>=2.1.0" -y
conda install -c anaconda graphviz -y

3行目はREADMEと違いますが、magattacaさんのブログ によるとGraphvizのインストール時には-c anacondaとするほうが良いそうなのでそうしました。（READMEでも言及されています）

AiZynthFinderのインストールはGitHubから最新版のソースコードを落とす方法を採用します。モデルをダウンロードするときに結局ソースコードを落とすことになるので先にコードをクローンしたほうが手間が省けます。

git clone https://github.com/MolecularAI/aizynthfinder.git
cd aizynthfinder
python -m pip install -e .

AiZynthFinderのインストールはこれで完了です。

動作確認

AiZynthFinderのインストールが終わると、aizynthcliとaizynthappの2つのコマンドが実行できるようになっています。動作確認を行うためにはまず学習済みモデルをダウンロードする必要があります。モデルをダウンロードするためのスクリプトがソースコードの aizynthfinder/tools/download_public_data.py に入っています。

ハマりポイント2

このスクリプトを実行するときにはダウンロード先のディレクトリを指定するのですが、ディレクトリを先に作っておく必要があります。地味なポイントですが、エラーメッセージが分かりにくいので混乱しました。次のようにすればいいです。

mkdir model
python aizynthfinder/tools/download_public_data.py model

動作確認はaizynthappを使って行います。

aizynthapp --config model/config.yml

を実行するとブラウザが立ち上がりJupyter Notebookが開きます。最初のセルを実行すると次のような画面が表示されます。

SMILESと書いてあるテキストボックスに適当なSMILESを入れると逆合成が実行されます。試しに片頭痛治療剤として使われる Lomerizine の逆合成を実行します。

WikipediaによるとロメリジンのSMILESは
COC1=C(C(=C(C=C1)CN2CCN(CC2)C(C3=CC=C(C=C3)F)C4=CC=C(C=C4)F)OC)OC
です。

SMILESを入力したらRun Searchを押して逆合成を実行します。経路探索が完了するまでしばらくかかります。

経路探索が終わったら、Show Reactionsを押して結果を表示します。このような結果が得られます。

Option 4で提示された以下のルートは創薬化学―有機合成からのアプローチで紹介されているものとほぼ同じでした。（参考書にはこのルートの最後の反応だけが記載されていました）

正しい経路を予測できたと言ってもよいのかもしれませんが、AiZynthFinderは特許データベースを利用して機械学習を行っているのでこのルートが学習データに含まれている可能性があります。実際にGoogle PatentでLomerizineと検索したところ、EP0159566A1に該当する反応の記述がありました。この反応がデータベースに含まれている可能性が高そうです。

データベースに含まれていないが、合成経路の正解が分かっている化合物を使わないと正しく評価できないと思うのですが、どなたかそのような化合物の探し方をご存知ないでしょうか……

逆合成ちゃんを支える技術

AiZynthFinderの紹介が終わったので、逆合成ちゃんの要素技術を紹介していきます

SMILES取得: Tweepy

逆合成ちゃんへのリプライは、Twitter APIの statuses/filter を使って、@retrosynthchan という文字列を含むツイートをほぼリアルタイムで取得しています。実装にはPythonでTwitterを扱うライブラリのTweepyを使っています。ライブラリの使い方の詳細はStreaming With Tweepy に書いてあります。

逆合成ちゃんへのリプライは@retrosynthchan SMILESという文字列になっていることを期待しています。そのため、ツイート全文を空白スペースで分割した2番目の文字列がSMILESになっているとします。

入力されたSMILESに対してRDKitのMolFromSmiles()を呼び出してNoneが返ってきたら、SMILESとして不正な文字列が与えられたことになるので「SMILESを認識できませんでした」とリプライして終了します。SMILESが正しい場合には、RDKitで生成した画像を一旦リプライします。これは逆合成の実行には時間がかかるため、そのツイートに反応したことを早めにユーザーに伝えるためです。エラー処理などの詳細を省くと次のようなコードになります。

class MyStreamListener(tweepy.StreamListener):
    def on_status(self, status):
        smiles = status.text.split(' ')[1]
        mol = Chem.MolFromSmiles(smiles)
        if mol is None:
            api.update_status(
              status='@{} SMILESを認識できませんでした'.format(
	          status.author.screen_name),
              in_reply_to_status_id=status.id)
            return
        else:
	    # 分子画像を生成してリプライ
            Draw.MolToFile(mol, 'mol.png')
            res = api.media_upload('mol.png'
            api.update_status(
              status='@{} {}の逆合成を実行します。しばらくお待ちください。'.format(
                  status.author.screen_name, smiles),
              media_ids=[res.media_id],
              in_reply_to_status_id=status.id)

            # 逆合成の実行
	    run_retrosynthesis(smiles)

逆合成の実行: AiZynthFinder

SMILESを取得したらAiZynthFinderを実行します。ここではaizynthcliを利用します。正しくは Python interface を使うべきだと思いますが、それは今後の課題とします。

aizynthcliの詳細はドキュメントの Command-line interface に書かれています。

smiles.txt に逆合成を実行したいSMILESを書き込んで

aizynthcli --config config_local.yml --smiles smiles.txt

を実行するとoutput.hdf5というファイルが生成されます。このファイルに対して以下のコードを実行すると、合成経路の画像が生成されます。

import pandas as pd
from aizynthfinder.analysis import ReactionTree

data = pd.read_hdf("output.hdf5", "table")
all_trees = data.trees.values
trees_for_first_target = all_trees[0]

for itree, tree in enumerate(trees_for_first_target):
    imagefile = f"route{itree:03d}.png"
    ReactionTree.from_dict(tree).to_image().save(imagefile)

画像の結合: ImageMagick

Twitterには4つまでしか画像を貼れないので、合成経路の画像を3つずつ結合した画像ファイルを作ります。3つずつにしたのはそれ以上つなげると見にくかったからという理由です。

画像の結合にはImageMagickを使いました。ImageMagickは脆弱性が多いため、使わない選択をしている会社もありますが、今回処理するのは自分たちで生成した画像なのであまり気にしなくてもいいでしょう。

あまり深いことを考えずに subprocess.run() でImageMagickのconvertを実行しています。3つずつ分ける処理には More Itertools の chunked を使っています。

images = glob.glob('route*.png')
for i, imgs in enumerate(more_itertools.chunked(sorted(images), 3)):
    files = ' '.join(imgs)
    proc = subprocess.run("convert -append {} result{}.png".format(files, i),
                          shell=True, text=True)

結果の返信: Tweepy

こうして得られた逆合成の結果をTweepyでリプライします。

filenames = glob.glob('result*.png')
media_ids = []
for filename in sorted(filenames)[:4]:
    res = api.media_upload(filename)
    media_ids.append(res.media_id)
api.update_status(status='@{} {}の逆合成結果です'.format(
    status.author.screen_name, smiles),
    media_ids=media_ids, in_reply_to_status_id=status.id)

このようにして、入力されたSMILESに対して逆合成解析を行った結果を最大12個リプライするbotを作ることができました。

プロセスの監視

TweepyでStreaming APIを使うと途中で落ちることがよくあります。そこで、プロセスが落ちていないか監視し、落ちていたら再起動するシェルスクリプトをcronで1分ごとに走らせています。

#!/bin/bash

isAlive=`ps ax | grep "retrosynthesis.py" | grep -v grep | wc -l`
if [ $isAlive = 0 ]; then
nohup python -u /home/user/retrosynthesis.py >> /home/user/log.txt 2>&1 &
fi

systemctlとかを使うのが真っ当なやり方なのかもしれませんが、分かりやすいのでいつもこうしています……

conda環境でcronジョブを実行する方法はcron job to run under conda virtual environmentを参考にしました。

今後の課題

ジョブキューの導入

今のところ逆合成の同時実行数を制御する仕組みがないので、同時にリクエストが飛んでくると何も考えずに同時実行し落ちる場合があります。CeleryやRQなどを使ってジョブを管理できるようにすることを考えています。

ストックの設定

デフォルトの設定ではZINCデータベースにある化合物に分解されるまで逆合成が実行されます。ZINCデータベースにはかなり複雑な化合物も収録されているため、かなり早い段階で逆合成が終わってしまう例がよく見られました。もう少し小さな化合物になるまで逆合成が進むように適切なストックを設定する必要がありそうです。日本で入手しやすい試薬の一覧みたいなデータベースをご存知の方は教えていただけると嬉しいです。

サーバー代

AiZynthFinderの処理はかなり重いので、ConoHa VPS の上から2番目のプランで実行しています。ConoHa VPSはAWS Lightsailよりも値段に対してCPU性能が良く、さくらのVPSと違って最低利用期間がないので採用しました。
月額3,280円はちょっと高いので、今のところAdvent Calendar期間が終わったら停止する予定でいます。もし需要があるようならGCPなどのより安いサーバーに移行し、Web版を作って広告を設置するなどしてなんとかペイできるようにする方法を考えたいです。

おわりに

今回はAiZynthFynderの使いかたを中心に技術的な面を紹介しましたが、次回は今年読んだ一番好きな論文2020 Advent Calendar 2020でAiZynthFinderの論文を読んで理論的な面について解説する予定です。皆様の投稿した分子の逆合成経路についても触れてみたいので是非遊んでみてください。

Discussion