似ていても　意外と小さい　タニモト係数

## はじめに
**化合物の類似度を評価する際に用いられる「タニモト係数」について解説していきます。**
本投稿では以下について述べていきます。

- タニモト係数の定義
- タニモト係数とタニモト距離の関係
- タニモト係数の算出方法
- フィンガープリントとタニモト係数の関係

※マテリアルズインフォマティクス関係の内容を他にも投稿していますので、よろしければ[こちら](https://qiita.com/oki_kosuke/items/34b1071224953776c211#%E6%8A%95%E7%A8%BF%E4%B8%80%E8%A6%A7)の一覧から他の投稿も見て頂けますと幸いです。

## タニモト係数とは
**化合物間の類似度を評価する指標です。**類似した化合物は似た性質を持つとされていますので、ある化合物が他の化合物とどれぐらい類似しているかを知ることは材料探索において有用であると考えられます。タニモト係数は化合物の分子構造をフィンガープリントに直し、フィンガープリント同士の類似度を計算した指標です。つまり、後述しますが**用いるフィンガープリントによってタニモト係数の値は変化**します。

化合物Aと化合物Bのタニモト係数は、**化合物AとBに共通する部分構造の数を化合物AとBに含まれる部分構造の総数で割った値です。**式で表すと次のようになります。n(A∩B)は共通する部分構造の数、n(A∪B)は部分構造の総数を表す自然数とします。


```math
c_{Tanimoto}(A,B) = \frac{n(A∩B)}{n(A∪B)}
```


具体例で考えてみます。化合物Aと化合物Bが以下のようなフィンガープリント（fp）を持っているとします。なお、フィンガープリントでは部分構造が存在するか否かを0/1で表します。fp1は化合物Aには含まれるが、化合物Bには含まれないことを表しています。

![fp.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/3812b83a-4709-a54b-89a8-2fdf63e41809.png)

この場合
n(A∩B)：4（fp2, 3, 5, 10)
n(A∪B)：7（fp1, 2, 3, 5, 6, 8, 10）
なので、タニモト係数は4/7=0.57となります。**タニモト係数が1に近いほど化合物間の類似度が高い**と言えます。

### タニモト距離とは
1からタニモト係数を引いた値です。式で表すと次のようになります。

```math
d_{Tanimoto}(A,B) = 1- c_{Tanimoto}(A,B)
```

### タニモト係数vsタニモト距離

- タニモト係数は類似度の指標
- タニモト距離は距離指標

表にまとめると以下の関係になります。

![tanimoto.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/e56a45c8-0f29-9f04-c94d-8e4975bca8aa.png)


## タニモト係数の算出

### 環境
- windows10
- conda 4.9.2
- python 3.7.1
- rdkit 2020.03.2.0

### 具体的なタニモト係数の算出方法
[こちら](https://qiita.com/oki_kosuke/items/51e461ec226eca9454bc)の記事で最大共通部分骨格（MCS）を求めた3分子（mol1, mol2, mol3）を対象にタニモト係数を算出していきます。この投稿からご覧の方のためにも3分子のデータの取得方法を以下に記載します。なお、データは[RDKitのGitHubページ](https://github.com/rdkit/rdkit/tree/master/Docs/Book/data)から取得した"cdk2.sdf"のデータを読み込んでいます。

```python
from rdkit import Chem

# データの読み込み
suppl = Chem.SDMolSupplier('data/cdk2.sdf')
ms = [x for x in suppl if x is not None]

# msのデータを共通骨格pを持つ分子に絞り、pでアライメントしておく
p = Chem.MolFromSmiles('[nH]1cnc2cncnc21')
subms = [x for x in ms if x.HasSubstructMatch(p)]
from rdkit.Chem import AllChem, Draw
AllChem.Compute2DCoords(p)
for m in subms: AllChem.GenerateDepictionMatching2DStructure(m,p)
img=Draw.MolsToGridImage(subms,molsPerRow=4,subImgSize=(200,200),legends=[x.GetProp("_Name") for x in subms])

# submsの5,6,7番目をリストとして保存
from rdkit.Chem import rdFMCS
mol1 = subms[5]
mol2 = subms[6]
mol3 = subms[7]
mols = [mol1, mol2, mol3]
Draw.MolsToGridImage(mols,molsPerRow=len(mols),subImgSize=(300,300))
```

![subms_wo_MSC.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/251b3661-f432-2ebe-d565-6e1b8fc7a561.png)

続いてフィンガープリントを計算し、タニモト係数を求めていきます。今回はmol1（左）を基準にして、mol1に対するタニモト係数を算出していきます。

```python
from rdkit import DataStructs

# フィンガープリントを定義する（radius=2, bitInfo=2048を指定）
morgan_fp = [AllChem.GetMorganFingerprintAsBitVect(mol, 2, 2048) for mol in mols]
#molsの0番目のmorgan_fpを基準にして、各タニモト係数(最小0～最大1)を計算する
tanimoto = DataStructs.BulkTanimotoSimilarity(morgan_fp[0], morgan_fp)
Draw.MolsToGridImage(mols, molsPerRow=3, subImgSize=(300,300), legends=['Tanimoto: {:.2f}'.format(i) for i in tanimoto])
```
<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/2de2ea04-a6ae-f040-eb99-9ae4b32aa59d.png" width="100%">

タニモト係数は順に**1.00, 0.52, 0.38**という結果でした。構造を見た感じではmol2（真ん中）、mol3（右）ともにmol1（左）とかなり類似していますが、タニモト係数は結構小さいなという印象を受けました。なお、morgan_fpはrdkit.DataStructs.cDataStructs.ExplicitBitVectがlistとして保存されています。


## 類似構造間のタニモト係数
構造の類似度とタニモト係数についてもう少し検討していきます。例として、アセトアミノフェンを基準に少しずつ構造を変更した分子のタニモト係数を求めてみます。タニモト係数の基準化合物はmol1のアセトアミノフェンとします。

```python
# アセトアミノフェンを定義
mol1 = Chem.MolFromSmiles('CC(=O)Nc1ccc(O)cc1')

# アセトアミノフェンから少しずつ構造を変化させた化合物を定義
# mol1の側鎖にメチル基を追加
mol2 = Chem.MolFromSmiles('CCC(=O)Nc1ccc(O)cc1')
# mol2の側鎖にメチル基を追加
mol3 = Chem.MolFromSmiles('CC(C)C(=O)Nc1ccc(O)cc1')
# mol1のC=OをC-OHに変更
mol4 = Chem.MolFromSmiles('CC(O)Nc1ccc(O)cc1')
# mol1のN-HのHをメチル基に置換
mol5 = Chem.MolFromSmiles('CC(=O)N(C)c1ccc(O)cc1')

# フィンガープリントを計算
list = [mol1, mol2, mol3, mol4, mol5]
morgan_fp = [AllChem.GetMorganFingerprintAsBitVect(mol, 2, 2048) for mol in list]
# molsの0番目のmorgan_fpを基準にして、各タニモト係数(最小0～最大1)を計算する
tanimoto = DataStructs.BulkTanimotoSimilarity(morgan_fp[0], morgan_fp)
Draw.MolsToGridImage(list, molsPerRow=len(list), subImgSize=(300,300), legends=['Tanimoto: {:.2f}'.format(i) for i in tanimoto])
```

<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/6f09473e-77ea-ca74-7fa2-67602bebb706.png" width="100%">

タニモト係数は順に**1.00, 0.59, 0.59, 0.43, 0.41**という結果でした。メチル基を側鎖に1つ追加することで、タニモト係数が1→0.59に変化することが分かりました。（感覚的には0.59よりもっと大きな値になるイメージでした）

### フィンガープリントを変えた際のタニモト係数
タニモト係数は定義からも分かるように、用いるフィンガープリントによってその値が変わります。今回はMorganフィンガープリントで半径（radius）を変えた際にタニモト係数がどのような値を取るのかについて検討していきます。なお、上記のコードではradius=2で計算したので、radius=3と1の場合を検討していきます。

##### radius=3の場合

```python
# radius=3としてフィンガープリントを定義
morgan_fp = [AllChem.GetMorganFingerprintAsBitVect(mol, 3, 2048) for mol in list]
#molsの0番目のmorgan_fpを基準にして、各タニモト係数(最小0～最大1)を計算する
tanimoto = DataStructs.BulkTanimotoSimilarity(morgan_fp[0], morgan_fp)
Draw.MolsToGridImage(list, molsPerRow=len(list), subImgSize=(300,300), legends=['Tanimoto: {:.2f}'.format(i) for i in tanimoto])
```

<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/0d5ea9e1-336e-b276-b988-db460fd88332.png" width="100%">

タニモト係数は順に**1.00, 0.54, 0.59, 0.39, 0.34**という結果でした。

##### radius=1の場合

```python
# radius=1としてフィンガープリントを定義
morgan_fp = [AllChem.GetMorganFingerprintAsBitVect(mol, 1, 2048) for mol in list]
#molsの0番目のmorgan_fpを基準にして、各タニモト係数(最小0～最大1)を計算する
tanimoto = DataStructs.BulkTanimotoSimilarity(morgan_fp[0], morgan_fp)
Draw.MolsToGridImage(list, molsPerRow=len(list), subImgSize=(300,300), legends=['Tanimoto: {:.2f}'.format(i) for i in tanimoto])
```

<img src="https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/83a0a9c0-dda4-2d49-7080-c2f6d91f9e46.png" width="100%">

タニモト係数は順に**1.00, 0.67, 0.67, 0.47, 0.53**という結果でした。

結果は以下の表の通りです。radiusを小さくするとより狭い範囲で構造を比較することになるので、タニモト係数が大きくなる傾向にあることは感覚的には理解できます。mol3はradius=2と3でタニモト係数が同じである点、mol4とmol5はradiusによってタニモト係数の大小関係が逆転している点は興味深いです。フィンガープリントとタニモト係数の関係については、また深堀していきたいと考えています。
![tanimoto2.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1057114/185c306f-7063-ab79-ada7-eef5a1190f42.png)


## まとめ

**化合物の類似度を評価する際に用いられる「タニモト係数」について解説いたしました。**
結論としては以下の通りです。

- タニモト係数はフィンガープリントを基に化合物の類似度を評価する指標である
- タニモト係数は類似度、タニモト距離は距離の指標であり相反する関係にある
- タニモト係数はrdkitのDataStructs.BulkTanimotoSimilarityで算出する
- 用いるフィンガープリントによってタニモト係数は変化する

### 参考
- [化学の新しいカタチ](https://future-chem.com/rdkit-fingerprint/)
- [あなたのおっしゃる「タニモト係数」って、どれのことかしら？](https://qiita.com/maskot1977/items/fac2931016aa7dbcbd6f)
- [実践マテリアルズインフォマティクス](https://www.amazon.co.jp/%E5%AE%9F%E8%B7%B5-%E3%83%9E%E3%83%86%E3%83%AA%E3%82%A2%E3%83%AB%E3%82%BA%E3%82%A4%E3%83%B3%E3%83%95%E3%82%A9%E3%83%9E%E3%83%86%E3%82%A3%E3%82%AF%E3%82%B9-%E8%88%B9%E6%B4%A5-%E5%85%AC%E4%BA%BA/dp/4764906155)
- [RDKitの部分構造についての理解を深める](https://qiita.com/oki_kosuke/items/51e461ec226eca9454bc)