指定した回数だけ化学反応を実施する関数の作成

今回は指定した回数だけ化学反応を実施して生成物を列挙する自作関数の作成を紹介します。
[以前の記事][9]で似たような繰り返し反応について紹介しましたが、今回は繰り返し回数に終点がなく
指定した回数だけ延々と続けられるような反応に向いた自作関数を紹介します。

---
###これまでの記事のシリーズ

[RDKit入門①：まずは環境構築][1]
[RDKit入門②：1分子の読み込みと部分構造検索][2]
[RDKit入門③：複数分子の読み込み (前編)と分子群の描画][3]
[RDKit入門④：複数分子の読み込み (後編)とデータフレームの加工][4]
[RDKit入門⑤：データフレーム内の分子群に対する部分構造検索][5]
[RDKit入門⑥：Morganフィンガープリントの作成とそれを用いたタニモト係数の計算による分子類似性評価][6]
[RDKit入門⑦：反応式の取り扱い(前編)][7]
[RDKitで化学反応を扱った際の生成物群を整理し図示する方法][8]
[RDKit入門⑧：反応式の取り扱い(後編)][9]
[RDKitで化学反応を繰り返しつつ分子の重複数をカウント・記載する][10]
[RDKit入門⑨：データフレーム内の分子群に対する化学反応の適用とExcelファイルへの出力][11]
[RDKit入門⑩：データフレーム内の分子群に対する完全一致検索][12]

[1]:https://qiita.com/Ru-PNP/items/ed5a8530611ef0d712b7
[2]:https://qiita.com/Ru-PNP/items/58d733ef327677d273fd
[3]:https://qiita.com/Ru-PNP/items/1cea34b864e993e4e1e1
[4]:https://qiita.com/Ru-PNP/items/62a108d92c1ddd522e58
[5]:https://qiita.com/Ru-PNP/items/b469ac1b907909a9d821
[6]:https://qiita.com/Ru-PNP/items/79873f0e3eab444f5bfe
[7]:https://qiita.com/Ru-PNP/items/49a080a015414553bbad
[8]:https://qiita.com/Ru-PNP/items/baebf99e269bd7fb2fb5
[9]:https://qiita.com/Ru-PNP/items/e0c1cccf4a7b0cd3f412
[10]:https://qiita.com/Ru-PNP/items/25d40044f3985a406982
[11]:https://qiita.com/Ru-PNP/items/897541250b69733f6f1e
[12]:https://qiita.com/Ru-PNP/items/6185eb2fe873e9d9a893

###環境
- Windows 10
- python 3.8.8
- RDKit 2021.03.1

---
###下準備
```python
from rdkit import Chem
from rdkit.Chem import AllChem, Draw
from itertools import chain
#パッケージの読み込み

mol_1 = Chem.MolFromSmiles("CO")
mol_2 = Chem.MolFromSmiles("NCCOCC(=O)Nc1cccc2c1C(=O)N(C1CCC(=O)NC1=O)C2=O")
rxn_1 = AllChem.ReactionFromSmarts("[Ch:1]>>[C:1]-C")
rxn_2 = AllChem.ReactionFromSmarts("[Nh2]CC[O:1][C:2]>>[NH2]CC[O:1]CCO[C:2]")
#今回の題材となる基質と反応をそれぞれ定義
```
今回はmethanol (mol_1)を出発にアルキル鎖を伸長する反応(rxn_1)と
Pomalidomide-PEG1-NH2 (mol_2)を出発にPEG鎖を伸長する反応(rxn_2)を題材とします。

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/e9b3c550-d81d-6d0a-9f4e-7cad4efed57f.png)

---
###指定した回数だけ化学反応を繰り返し生成物を網羅する自作関数
原料(mol)・反応実施数(n)・反応(rxn)の3つを指定する関数を作成しました。
結果はリストとして出力されます。
なおリストには原料(mol)も出力されるように指定しましたが、必要に応じて書き換えてください。

```python
def repeat_reactions(mol, n, rxn):
    """
    Parameters
    ----------
    mol : Mol object of rdkit.Chem.rdchem module
        繰り返し反応の出発物質を入力する
    n : int
        反応を繰り返す数を入力する
    rxn : ChemicalReacion object of rdkit.Chem.rdChemReactions module
        繰り返したい反応を入力する

    Returns
    -------
    all_products : list
        反応をn回まで繰り返す際に登場したmolオブジェクトが全て格納されたリストが出力される
    """

    #ループ前の処理
    rxn.Initialize()    #念のため反応を初期化しておく。
    mol.UpdatePropertyCache(strict = False) #念のため原料の価数等を更新しておく。
    all_products = [mol]    #原料の格納
    products_unique = [mol] #一周目の反応の原料を格納
    i = 0
    
    #ループ中の処理
    while i < n:
        products = []   #空のリストを用意
        for product in products_unique:
            if rxn.IsMoleculeReactant(product):
               products.extend(list(chain.from_iterable(rxn.RunReactants([product]))))  
            #i-1サイクル目の反応の生成物がiサイクル目の反応の原料となるかどうかを確認。
            #原料となる場合は反応を実施し、空のリストproductsに格納。                
        products_unique = [Chem.MolFromSmiles(smile) for smile in {Chem.MolToSmiles(x) for x in products}]
        #productsの中の分子群中の重複を削除し、products_uniqueという名のリストに格納。
        for product in products_unique:
            product.UpdatePropertyCache(strict = False)
        #iサイクル目の反応の生成物の価数等を更新しておく
        all_products.extend(products_unique)    #all_productsにproducts_uniqueの中身を追加
        i += 1  #次のサイクルへ
    
    #ループが終了した際の処理
    else:
        return all_products
```
mol_1を出発原料としてrxn_1を4回適用した例を示します。

```python
result_1 = repeat_reactions(mol_1, 4, rxn_1)
Draw.MolsToGridImage(result_1, molsPerRow = 4)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/9f626002-1f46-665a-3a1e-800aaffd40b6.png)
同様にmol_2を出発原料としてrxn_2を5回適用した例を示します。

```python
result_2 = repeat_reactions(mol_2, 5, rxn_2)
Draw.MolsToGridImage(result_2, subImgSize = (300,300))
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/9562f886-2dee-b241-8970-f7513ea7af51.png)

出力結果は適当な情報と共にExcelファイルに出力してもよいかもしれません。
（Excelファイルへの出力については[前回の記事][12]を参照）

```python
df = pd.DataFrame(result_2, columns = ['ROMol'])
#result_2からのデータフレームの作成

df['分子式'] = [rdMolDescriptors.CalcMolFormula(x) for x in result_2]
df['分子量'] = [round(Descriptors.MolWt(x), 2) for x in result_2]
df['MolLogP'] = [Descriptors.MolLogP(x) for x in result_2]
df['水素結合ドナー数'] = [rdMolDescriptors.CalcNumHBD(x) for x in result_2]
df['水素結合アクセプター数'] = [rdMolDescriptors.CalcNumHBA(x) for x in result_2]
df['SMILES'] = [Chem.MolToSmiles(x) for x in result_2]
#記述子やSMILES構造式の追加

PandasTools.SaveXlsxFromFrame(df, "result_2.xlsx", molCol = 'ROMol')
#Excelファイルへの出力
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/14fe348f-0de5-fa69-4ab6-07b32f6901e5.png)


---
###最後のサイクルに登場する化合物（群）だけを取り出す自作関数
先に紹介した関数では原料および1～n回目の反応で登場する化合物を網羅しましたが、
n回目の結果だけ欲しい場合についても作成しました。
結果は分子数が1つの場合はmolオブジェクト、複数の場合はリストで出力されるようにしました。

```python
def repeat_reactions_nth_cycle(mol, n, rxn):
    """
    Parameters
    ----------
    mol : Mol object of rdkit.Chem.rdchem module
        繰り返し反応の出発物質を入力する
    n : int
        反応を繰り返す数を入力する
    rxn : ChemicalReacion object of rdkit.Chem.rdChemReactions module
        繰り返したい反応を入力する

    Returns
    -------
    products_unique : list
        nサイクル目の反応で登場するmolオブジェクトが格納されたリストが出力される。
        ただしproducts_uniqueの中身が1つの場合は中身products_unique[0]を出力する。
    """
    #ループ前の処理
    rxn.Initialize()    #念のため反応を初期化しておく。
    mol.UpdatePropertyCache(strict = False) #念のため原料の価数等を更新しておく。
    all_products = [mol]    #原料の格納
    products_unique = [mol] #一周目の反応の原料を格納
    i = 0
    
    #ループ中の処理
    while i < n:
        products = []   #空のリストを用意
        for product in products_unique:
            if rxn.IsMoleculeReactant(product):
               products.extend(list(chain.from_iterable(rxn.RunReactants([product]))))  
            #i-1サイクル目の反応の生成物がiサイクル目の反応の原料となるかどうかを確認。
            #原料となる場合は反応を実施し、空のリストproductsに格納。                
        products_unique = [Chem.MolFromSmiles(smile) for smile in {Chem.MolToSmiles(x) for x in products}]
        #productsの中の分子群中の重複を削除し、products_uniqueという名のリストに格納。
        for product in products_unique:
            product.UpdatePropertyCache(strict = False)
        #iサイクル目の反応の生成物の価数等を更新しておく
        all_products.extend(products_unique)    #all_productsにproducts_uniqueの中身を追加
        i += 1  #次のサイクルへ
    
    #ループが終了した際の処理
    else:
        if len(products_unique) == 1:
            return products_unique[0]
        else:
            return products_unique
```
mol_1を出発原料としてrxn_1を4回適用した例を示します。
この例では結果がリストとして出力されます。

```python
result_3 = repeat_reactions_nth_cycle(mol_1, 4, rxn_1)
Draw.MolsToGridImage(result_3, molsPerRow = 4)
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/c1a95943-9678-7cb8-221f-d11736947436.png)

mol_2を出発原料としてrxn_2を5回適用した例を示します。
この例では結果はmolオブジェクトとして出力されます。

```python
result_4 = repeat_reactions_nth_cycle(mol_2, 5, rxn_2)
result_4
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1144410/3a81d985-546b-4ac2-793b-408d8add06eb.png)

---
###おわりに
今回は指定した回数だけ化学反応を実施する自作関数を作成いたしました。
扱う反応を上手に選べば多様な化合物ライブラリーを自動的に作成でき、
いちいちSMILES構造式を打ち込む手間を省略することができます。
構造活性相関表の作成等、使用できる場面があると思いますでご活用ください。

###参考
- [RDKitにおける記述子の扱い方をリピンスキーの法則を通して学ぶ][13]

[13]:https://future-chem.com/rdkit-descriptors/