【医薬品業界】MF公示からデータ分析した件

はじめに

この記事は医薬品業界の中でも後発医薬品 (ジェネリック)に関する記事です。
データ分析のアプローチとして、こんな方法もあるんだな程度で見ていただければと思います。

目次

• 背景
• 目的
• 方法
• 結果・考察
• 所感
• 今後の予定

背景

私はジェネリック医薬品企業に勤める者です。恥ずかしながら、競合他社がどこなのか把握していませんでした。
一方、PMDA (医薬品機構)はMF公示資料を月に2度更新しています (リンク)。MF (Master File)とは、医薬品の原薬 (有効成分)を製造する企業がPMDAに登録する製造方法のノウハウのことです。このMFを医薬品の製剤を製造する企業が参照することで、原薬製造企業はノウハウを開示することなく、医薬品が作ることができます。
つまり、このMF公示資料を分析すれば、競合他社を把握することができるのでは？と考えました。

目的

以下2つのタスクを考えました。
競合他社の動向調査：他社がどのくらいアクティブかを分析する
競合他社の技術力調査：他社がどのくらい合成難度の高い化合物に挑戦しているかを分析する

方法

競合他社の動向調査

1. これまで収集したMF公示資料 (excel)をpandasで一つのDataFrameにする。
2. 横軸に年代、縦軸に会社名、バブルの大きさをMF件数とするバブルチャートを作成する。

競合他社の技術力調査

1. 上記で1つにしたDataFrameを使用する。
2. 化合物名+casで検索してcas番号の列を挿入する。
3. cas番号からpubchemでsmilesを取得後、rdkitでMolObjectを入手する。
4. MolObjectからSynthetic Accessibility score (sa)を算出し、平均及び標準偏差を取得する。

実行環境は以下の通りです。

• python 3.7.11
• pandas 1.3.4
• numpy 1.20.3
• rdkit 2020.09.1
• matplotlib 3.4.2
• requests 2.27.1
• re 2.2.1
• pubchempy 1.0.4

結果・考察

競合他社の動向調査

まずはライブラリのインポートです。フォント指定をしないと、日本語が文字化けするため、plt.rcParams['font.family'] = 'MS Gothic'でフォント指定しています。

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from glob import glob

plt.rcParams['font.family'] = 'MS Gothic'

これまで収集したMF公示資料をMF_dataというフォルダに保存し、それを結合しました。メモリ使用量を削減するために、結合するごとに重複を除いていきました。

for i, file in enumerate(glob('MF_data/*.xlsx')):
    df = pd.read_excel(file)
    if i==0:
        df_con = df
        continue
    df_con = pd.concat([df_con, df])
    df_con = df_con.drop_duplicates()
df_con.head()

"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所
0	231MF10089	2019/04/11	2019/04/11	医薬品等原薬	セフメタゾール	Ａｍｉｃｏｇｅｎ（Ｃｈｉｎａ）Ｂｉｏｐｈａｒｍ　Ｃｏ．，Ｌｔｄ．	山東省済寧市詩仙路１６８８号	NaN	NaN
1	231MF10088	2019/04/11	2019/04/11	医薬品等原薬	アンブロキソール塩酸塩	ＡＬＣＨＹＭＡＲＳ　ＩＣＭ　ＳＭ　ＰＲＩＶＡＴＥ　ＬＩＭＩＴＥＤ．	“ＡＤＷＡＶＥ　ＴＯＷＥＲＳ”　Ｃ／４，ＴＨＩＲＤ　ＦＬＯＯＲ　ＮＥＷ　ＮＯ．１７（ＯＬＤ　...	NaN	NaN
2	231MF10087	2019/04/11	2019/04/11	医薬品等原薬	ソリスロマイシン	富士フイルムワコーケミカル株式会社	宮崎県宮崎市清武町木原３５５８番地	NaN	NaN
3	231MF10086	2019/04/11	2019/04/11	医薬品等原薬	ブリンゾラミド	Ｉｃｒｏｍ　Ｓ．ｐ．Ａ．	Ｖｉａ　Ｄｅｌｌｅ　Ａｒｔｉ　３３，２０８６３　Ｃｏｎｃｏｒｅｚｚｏ　（ＭＢ），　Ｉｔａｌｙ	NaN	NaN
4	227MF10236	2015/09/24	2019/04/11	医薬品等原薬	ゾルピデム酒石酸塩	Ｇｌｅｎｍａｒｋ　Ｌｉｆｅ　Ｓｃｉｅｎｃｅｓ　Ｌｉｍｉｔｅｄ	Ｇｌｅｎｍａｒｋ　Ｈｏｕｓｅ，ＨＤＯ　Ｃｏｒｐｏｒａｔｅ　Ｂｕｉｌｄｉｎｇ，Ｗｉｎｇ　Ａ，Ｂ...	NaN	NaN
"""

結合したDataFrameの詳細を確認します。

print(df_con.shape)
print(df_con['登録区分'].unique())
df_con.info()

"""
(6552, 9)
['医薬品等原薬' '添加剤' '医療機器原材料' 'その他']
<class 'pandas.core.frame.DataFrame'>
Int64Index: 6552 entries, 0 to 4306
Data columns (total 9 columns):
 #   Column       Non-Null Count  Dtype 
---  ------       --------------  ----- 
 0   MF登録番号       6552 non-null   object
 1   初回登録年月日      6552 non-null   object
 2   最新登録年月日      6552 non-null   object
 3   登録区分         6552 non-null   object
 4   登録品目名        6552 non-null   object
 5   登録者氏名        6552 non-null   object
 6   登録者住所        6552 non-null   object
 7   原薬等国内管理人の名称  33 non-null     object
 8   原薬等国内管理人の住所  33 non-null     object
dtypes: object(9)
"""

ここまでで分かったことを列挙します。

• 6552行、9列のデータである。
• MF登録番号は通し番号である。
• 初回登録年月日と最新登録年月日はよくわからなかったが、最新登録年月日が申請時の年月日である。
• 登録区分は、医薬品等原薬、添加剤、医療機器原材料、その他がある。
• 登録品目名は、原薬の化合物名が記されている。
• 登録者氏名は、登録した会社名である。
• 列名は「MF登録番号、初回登録年月日、最新登録年月日、登録区分 登録品目名、登録者氏名、登録者住所、原薬等国内管理人の名称、原薬等国内管理人の住所」
• 原薬等国内管理人の名称、原薬等国内管理人の住所は33行しか値がないので、情報としては価値がない

バブルチャートを作る準備をします。こちらを参考にして、行を会社名、列を年数、データを申請数としたDataFrameを作成しました。
初回登録年月日と最新登録年月日をdatetimeに変換し、そこから年だけを抽出してyearという列を生成しました。新たなDataFrameとして、行に登録者氏名、列にyearを持つdf_anaを作り、順次申請数を入力していきました。最後に、yearが年代順になっていなかったのでソートして下準備が完了しました。

df_con["初回登録年月日"] = pd.to_datetime(df_con["初回登録年月日"])
df_con["最新登録年月日"] = pd.to_datetime(df_con["最新登録年月日"])

df_con['year'] = df_con["最新登録年月日"].dt.year

df_ana = pd.DataFrame(index=df_con['登録者氏名'].unique(), columns=df_con['year'].unique())
for reg in df_ana.index:
    for year in df_ana.columns:
        a = (df_con['year']==year)&(df_con['登録者氏名']==reg)
        df_ana.loc[reg, year] = a.values.sum()

df_ana.columns = df_ana.columns.astype(int)
df_ana.reindex(columns = list(range(2005, 2023)))

df_ana.head()

	2005	2006	2007	2008	2009	2010	2011	2012	2013	2014	2015	2016	2017	2018	2019	2020	2021	2022
Ａｍｉｃｏｇｅｎ（Ｃｈｉｎａ）Ｂｉｏｐｈａｒｍ\u3000Ｃｏ．，Ｌｔｄ．	0	0	0	0	0	0	0	0	0	0	0	1	1	0	1	0	0	0
ＡＬＣＨＹＭＡＲＳ\u3000ＩＣＭ\u3000ＳＭ\u3000ＰＲＩＶＡＴＥ\u3000ＬＩＭＩＴＥＤ．	0	1	0	0	0	0	0	0	0	0	0	1	0	0	1	0	0	0
富士フイルムワコーケミカル株式会社	0	4	0	0	0	0	0	0	0	0	1	2	1	0	2	1	0	0
Ｉｃｒｏｍ\u3000Ｓ．ｐ．Ａ．	0	0	0	0	0	1	0	0	0	1	1	0	1	1	2	4	0	0
Ｇｌｅｎｍａｒｋ\u3000Ｌｉｆｅ\u3000Ｓｃｉｅｎｃｅｓ\u3000Ｌｉｍｉｔｅｄ	0	0	0	0	0	0	1	0	3	0	4	3	1	1	7	2	6	2

いよいよバブルチャートの作成です。全ての企業をプロットするのは大変なので、Top20を表示することにしました。
まず、申請数の合計をsumという列に代入し、それをもとにソートして、Top20を抽出しました。x軸に年数、y軸に会社名を記載したバブルチャートを作成しました。plt.scatter()のsが初期値ではバブルの大きさが小さかったので、50倍にしています。申請数をplt.text()でバブルの真ん中に記載しました。

df_ana['sum'] = df_ana.sum(axis=1)
df_ana = df_ana.sort_values('sum', ascending=False)

df_extract = df_ana.drop('sum', axis=1).head(20)

plt.figure(figsize=(10,10))
for n,y in enumerate(df_extract.columns):
    for m,x in enumerate(df_extract.index):
        num = int(df_extract[y][x])
        plt.scatter(y, x, s=num*50)
        if num!=0:
            plt.text(y, x, str(num), ha='center', va='center', c='white')
plt.xlim(df_extract.columns.min()-1, df_extract.columns.max()+1)
plt.xticks(np.arange(df_extract.columns.min()-1,
df_extract.columns.max()+1, 1), rotation=45)
plt.xlabel('year')
plt.title('Top 20 companies applying for MF')
plt.show()

さらに、登録区分について、円グラフを作成しました。こちらは確認程度だったので、グラフは示しませんが、Ｌｉｆｅ　Ｔｅｃｈｎｏｌｏｇｉｅｓ　Ｃｏｒｐｏｒａｔｉｏｎ以外は医薬品等原薬でした。当会社はその他の区分となっており、医薬品原薬の製造会社ではないことが分かりました。

fig, axes = plt.subplots(5,4, figsize=(15,15))
section = df_con['登録区分'].unique()
for ax, com, i in zip(axes.ravel(), df_extract.index, range(20)):
    df_life = df_con[df_con['登録者氏名']==com]
    cnt_list = []
    for item in section:
        try:
            cnt_list.append(df_life['登録区分'].value_counts()[item])
        except:
            cnt_list.append(0)
    ax.pie(cnt_list, labels=section, startangle=90, counterclock=False, autopct="%1.1f%%")
    if len(com)>=15:
        com = com[:15]
    ax.set_title("{}. ".format(i+1) + com)

Top20のランキングです。

top20_com = pd.DataFrame(df_extract.index)
top20_com
"""
0
0	Ｌｉｆｅ　Ｔｅｃｈｎｏｌｏｇｉｅｓ　Ｃｏｒｐｏｒａｔｉｏｎ
1	株式会社　三洋化学研究所
2	Ｓｉｃｏｒ　Ｓｏｃｉｅｔａ（｀）　Ｉｔａｌｉａｎａ　Ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　Ｓ．ｒ．ｌ．
3	アルプス薬品工業株式会社
4	Ｏｌｏｎ　Ｓ．ｐ．Ａ．
5	ダイト株式会社
6	Ｋｙｏｎｇｂｏ　Ｐｈａｒｍａｃｅｕｔｉｃａｌ　Ｃｏ．，Ｌｔｄ．
7	Ｍｙｌａｎ　Ｌａｂｏｒａｔｏｒｉｅｓ　Ｌｉｍｉｔｅｄ
8	住友化学株式会社
9	株式会社パーマケム・アジア
10	Ｃａｍｂｒｅｘ　Ｐｒｏｆａｒｍａｃｏ　Ｍｉｌａｎｏ　Ｓ．ｒ．ｌ．
11	Ｈｅｔｅｒｏ　Ｌａｂｓ　Ｌｉｍｉｔｅｄ
12	ＫＯＬＯＮ　ＬＩＦＥ　ＳＣＩＥＮＣＥ，ＩＮＣ．
13	ＳＡＮＯＦＩ　ＣＨＩＭＩＥ
14	Ａｕｒｏｂｉｎｄｏ　Ｐｈａｒｍａ　Ｌｉｍｉｔｅｄ
15	ＨＡＮＭＩ　ＦＩＮＥ　ＣＨＥＭＩＣＡＬ　ＣＯ．，ＬＴＤ．
16	相模化成工業株式会社
17	株式会社　大阪合成有機化学研究所
18	大和薬品工業株式会社
19	八代製薬株式会社
"""

考察

今回の調査で得た考察を以下に示します。

• MF申請数の合計でランキングしたところ、1位のLife technology (現在はthermo fisher scienticの一部)のMF内容はその他となっており、医薬品原薬を製造している企業ではなかった。
• 上位の半数は海外勢であり、企業の競争レベルは世界規模に及んでいることが読み取れた。
• 国内で勢力が衰えているのはパーマケムアジアのみであり、一方で海外勢の申請数は軒並み上昇傾向であった。
• 三洋化学研究所はジェネリック医薬品原薬の製造及び受託製造の企業。
• アルプス薬品工業は医薬品原薬だけでなく、食品添加物や化粧品も手掛ける企業。
• ダイトは医薬品原薬だけでなく、製剤事業も展開しており、製販の立場を理解できるアドバンテージがあると考えられる。

競合他社の技術力調査

まずはライブラリのインポートです。

import requests
import re
import pandas as pd
import numpy as np
import pubchempy as pcp
from rdkit import Chem
import sascorer
import matplotlib.pyplot as plt

plt.rcParams['font.family'] = 'MS Gothic'

まずは、上記で作成したDataFrameからTop20の企業のみを抽出しました。1129行まで絞り込みました。

df_test = df_con.query('登録者氏名 in {}'.format(list(top20_com.iloc[:,0])))
df_test
"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所
26	219MF10293	2007/09/10	2019/03/29	医薬品等原薬	スルバクタムナトリウム（製造専用）	Ｋｙｏｎｇｂｏ　Ｐｈａｒｍａｃｅｕｔｉｃａｌ　Ｃｏ．，Ｌｔｄ．	１７４，Ｓｉｒｏｋ－ｒｏ，Ａｓａｎ－ｓｉ，Ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ＫＯＲＥＡ	NaN	NaN
28	218MF10710	2006/08/14	2019/03/29	医薬品等原薬	日本薬局方　スルバクタムナトリウム（製造専用）	Ｋｙｏｎｇｂｏ　Ｐｈａｒｍａｃｅｕｔｉｃａｌ　Ｃｏ．，Ｌｔｄ．	１７４，Ｓｉｒｏｋ－ｒｏ，Ａｓａｎ－ｓｉ，Ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ＫＯＲＥＡ	NaN	NaN
36	218MF10971	2006/12/01	2019/03/19	医薬品等原薬	ナロキソン塩酸塩水和物	ＳＡＮＯＦＩ　ＣＨＩＭＩＥ	８２　Ａｖｅｎｕｅ　Ｒａｓｐａｉｌ　９４２５０　Ｇｅｎｔｉｌｌｙ　Ｆｒａｎｃｅ	NaN	NaN
42	230MF10133	2018/09/14	2019/03/12	医薬品等原薬	メチルプレドニゾロンコハク酸エステルナトリウム	ＳＡＮＯＦＩ　ＣＨＩＭＩＥ	８２　Ａｖｅｎｕｅ　Ｒａｓｐａｉｌ　９４２５０　Ｇｅｎｔｉｌｌｙ　Ｆｒａｎｃｅ	NaN	NaN
52	220MF10128	2008/05/15	2019/03/05	医薬品等原薬	メチルプレドニゾロンコハク酸エステル	ＳＡＮＯＦＩ　ＣＨＩＭＩＥ	８２　Ａｖｅｎｕｅ　Ｒａｓｐａｉｌ　９４２５０　Ｇｅｎｔｉｌｌｙ　Ｆｒａｎｃｅ	NaN	NaN
...	...	...	...	...	...	...	...	...	...
2633	219MF10105	2007/03/23	2012/03/30	医薬品等原薬	イトラコナゾール	Ｍｙｌａｎ　Ｌａｂｏｒａｔｏｒｉｅｓ　Ｌｉｍｉｔｅｄ	Ｐｌｏｔ　Ｎｏ．５６４／Ａ／２２，　Ｒｏａｄ　Ｎｏ．９２，　Ｊｕｂｉｌｅｅ　Ｈｉｌｌｓ，　Ｈ...	NaN	NaN
2821	222MF10253	2010/10/28	2010/10/28	医薬品等原薬	オロパタジン塩酸塩	ＫＯＬＯＮ　ＬＩＦＥ　ＳＣＩＥＮＣＥ，ＩＮＣ．	１１０，ＭＡＧＯＫＤＯＮＧ－ＲＯ，ＧＡＮＧＳＥＯ－ＧＵ，ＳＥＯＵＬ，ＫＯＲＥＡ	NaN	NaN
3478	217MF10151	2005/05/31	2008/08/18	医薬品等原薬	マレイン酸セチプチリン「サンヨー」	株式会社　三洋化学研究所	大阪府堺市美原区多治井１４８－１	NaN	NaN
3793	219MF10159	2007/05/15	2007/05/15	医薬品等原薬	チクロピジン塩酸塩「テバ」	Ｓｉｃｏｒ　Ｓｏｃｉｅｔａ（｀）　Ｉｔａｌｉａｎａ　Ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　Ｓ．ｒ．ｌ．	Ｐｉａｚｚａｌｅ　Ｌｕｉｇｉ　Ｃａｄｏｒｎａ，　４　－　２０１２３　ＭＩＬＡＮＯ，　Ｉｔａｌｙ	NaN	NaN
3809	218MF10270	2006/02/17	2007/04/27	医薬品等原薬	日本薬局方　メキタジン	住友化学株式会社	東京都中央区日本橋二丁目７番１号	NaN	NaN
1129 rows × 9 columns
"""

効果があるかは不明ですが、文字を小文字にしてスペースを削除しました。

df_test['登録品目名'] = df_test['登録品目名'].str.lower()
df_test['登録者氏名'] = df_test['登録者氏名'].str.lower().replace('　', " ")
df_test['登録者住所'] = df_test['登録者住所'].str.lower().replace('　', " ")

第一の関門は、登録品目名に記載されている化合物名からcas番号を引っ張り出すステージです。rtn_casという関数を作成しました。登録品目名に記載されている化合物名を受け取り、「化合物名+cas」でgoogle検索を行い、ページ中にあるcas番号を正規表現によって抽出するものです。処理にかなり時間がかかりましたが、何とか実行することができました。化合物によっては、cas番号を抽出できませんでしたが、次に進みました。

def rtn_cas(name):
    name_searched = name + " cas"
    num = 1 # pages that you want to search
    url = "https://www.google.co.jp/search?hl=ja&num={}&q={}".format(num, name_searched)
    res = requests.get(url)
    if res.status_code != 200:
        return False
    txt = res.text
    ptn = r"\d+-\d\d-\d" # cas number
    result = re.findall(ptn, txt)
    if len(result)==0:
        return False
    return result[0]

df_test['cas'] = df_test['登録品目名'].map(rtn_cas)
df_test

"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所	cas
26	219MF10293	2007/09/10	2019/03/29	医薬品等原薬	スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	69388-84-7
28	218MF10710	2006/08/14	2019/03/29	医薬品等原薬	日本薬局方　スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	False
36	218MF10971	2006/12/01	2019/03/19	医薬品等原薬	ナロキソン塩酸塩水和物	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	51481-60-8
42	230MF10133	2018/09/14	2019/03/12	医薬品等原薬	メチルプレドニゾロンコハク酸エステルナトリウム	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	2375-03-3
52	220MF10128	2008/05/15	2019/03/05	医薬品等原薬	メチルプレドニゾロンコハク酸エステル	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	False
...	...	...	...	...	...	...	...	...	...	...
2633	219MF10105	2007/03/23	2012/03/30	医薬品等原薬	イトラコナゾール	ｍｙｌａｎ　ｌａｂｏｒａｔｏｒｉｅｓ　ｌｉｍｉｔｅｄ	ｐｌｏｔ　ｎｏ．５６４／ａ／２２，　ｒｏａｄ　ｎｏ．９２，　ｊｕｂｉｌｅｅ　ｈｉｌｌｓ，　ｈ...	NaN	NaN	False
2821	222MF10253	2010/10/28	2010/10/28	医薬品等原薬	オロパタジン塩酸塩	ｋｏｌｏｎ　ｌｉｆｅ　ｓｃｉｅｎｃｅ，ｉｎｃ．	１１０，ｍａｇｏｋｄｏｎｇ－ｒｏ，ｇａｎｇｓｅｏ－ｇｕ，ｓｅｏｕｌ，ｋｏｒｅａ	NaN	NaN	False
3478	217MF10151	2005/05/31	2008/08/18	医薬品等原薬	マレイン酸セチプチリン「サンヨー」	株式会社　三洋化学研究所	大阪府堺市美原区多治井１４８－１	NaN	NaN	False
3793	219MF10159	2007/05/15	2007/05/15	医薬品等原薬	チクロピジン塩酸塩「テバ」	ｓｉｃｏｒ　ｓｏｃｉｅｔａ（｀）　ｉｔａｌｉａｎａ　ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　ｓ．ｒ．ｌ．	ｐｉａｚｚａｌｅ　ｌｕｉｇｉ　ｃａｄｏｒｎａ，　４　－　２０１２３　ｍｉｌａｎｏ，　ｉｔａｌｙ	NaN	NaN	False
3809	218MF10270	2006/02/17	2007/04/27	医薬品等原薬	日本薬局方　メキタジン	住友化学株式会社	東京都中央区日本橋二丁目７番１号	NaN	NaN	False
1129 rows × 10 columns
"""

次の関門は、cas番号からsmilesを得るステージです。pubchempyを用いて、smilesを取得しました。pubchempyはPubChemから情報を取り出せるライブラリです。cas2smilesという関数を作成し、smilesという列に代入しました。

def cas2smiles(cas):
    if cas == False:
        return False
    result = pcp.get_properties(['CanonicalSMILES'], cas, 'name')
    if len(result)==0:
        return False
    smiles = result[0]['CanonicalSMILES']
    return smiles

df_test['smiles'] = df_test['cas'].map(cas2smiles)
df_test

"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所	cas	smiles
26	219MF10293	2007/09/10	2019/03/29	医薬品等原薬	スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	69388-84-7	CC1(C(N2C(S1(=O)=O)CC2=O)C(=O)[O-])C.[Na+]
28	218MF10710	2006/08/14	2019/03/29	医薬品等原薬	日本薬局方　スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	False	False
36	218MF10971	2006/12/01	2019/03/19	医薬品等原薬	ナロキソン塩酸塩水和物	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	51481-60-8	C=CCN1CCC23C4C(=O)CCC2(C1CC5=C3C(=C(C=C5)O)O4)...
42	230MF10133	2018/09/14	2019/03/12	医薬品等原薬	メチルプレドニゾロンコハク酸エステルナトリウム	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	2375-03-3	CC1CC2C3CCC(C3(CC(C2C4(C1=CC(=O)C=C4)C)O)C)(C(...
52	220MF10128	2008/05/15	2019/03/05	医薬品等原薬	メチルプレドニゾロンコハク酸エステル	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	False	False
...	...	...	...	...	...	...	...	...	...	...	...
2633	219MF10105	2007/03/23	2012/03/30	医薬品等原薬	イトラコナゾール	ｍｙｌａｎ　ｌａｂｏｒａｔｏｒｉｅｓ　ｌｉｍｉｔｅｄ	ｐｌｏｔ　ｎｏ．５６４／ａ／２２，　ｒｏａｄ　ｎｏ．９２，　ｊｕｂｉｌｅｅ　ｈｉｌｌｓ，　ｈ...	NaN	NaN	False	False
2821	222MF10253	2010/10/28	2010/10/28	医薬品等原薬	オロパタジン塩酸塩	ｋｏｌｏｎ　ｌｉｆｅ　ｓｃｉｅｎｃｅ，ｉｎｃ．	１１０，ｍａｇｏｋｄｏｎｇ－ｒｏ，ｇａｎｇｓｅｏ－ｇｕ，ｓｅｏｕｌ，ｋｏｒｅａ	NaN	NaN	False	False
3478	217MF10151	2005/05/31	2008/08/18	医薬品等原薬	マレイン酸セチプチリン「サンヨー」	株式会社　三洋化学研究所	大阪府堺市美原区多治井１４８－１	NaN	NaN	False	False
3793	219MF10159	2007/05/15	2007/05/15	医薬品等原薬	チクロピジン塩酸塩「テバ」	ｓｉｃｏｒ　ｓｏｃｉｅｔａ（｀）　ｉｔａｌｉａｎａ　ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　ｓ．ｒ．ｌ．	ｐｉａｚｚａｌｅ　ｌｕｉｇｉ　ｃａｄｏｒｎａ，　４　－　２０１２３　ｍｉｌａｎｏ，　ｉｔａｌｙ	NaN	NaN	False	False
3809	218MF10270	2006/02/17	2007/04/27	医薬品等原薬	日本薬局方　メキタジン	住友化学株式会社	東京都中央区日本橋二丁目７番１号	NaN	NaN	False	False
1129 rows × 11 columns
"""

ここまででどのくらい情報が欠損したのかチェックしました。1129件中408件を抽出できていました。かなり欠損していますが、仕方がないのでこのまま進めました。

(df_test['smiles']==False).value_counts()
"""
True     721
False    408
Name: smiles, dtype: int64
"""

そして、smilesからMolObjectを返す関数smiles2molを作成して、MolObjectを取得しました。

def smiles2mol(smiles):
    if smiles==False:
        return False
    mol = Chem.MolFromSmiles(smiles)
    return mol

df_test['ROMol'] = df_test['smiles'].apply(smiles2mol)
df_test
"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所	cas	smiles	ROMol
26	219MF10293	2007/09/10	2019/03/29	医薬品等原薬	スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	69388-84-7	CC1(C(N2C(S1(=O)=O)CC2=O)C(=O)[O-])C.[Na+]	<img data-content="rdkit/molecule" src="data:i...
28	218MF10710	2006/08/14	2019/03/29	医薬品等原薬	日本薬局方　スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	False	False	False
36	218MF10971	2006/12/01	2019/03/19	医薬品等原薬	ナロキソン塩酸塩水和物	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	51481-60-8	C=CCN1CCC23C4C(=O)CCC2(C1CC5=C3C(=C(C=C5)O)O4)...	<img data-content="rdkit/molecule" src="data:i...
42	230MF10133	2018/09/14	2019/03/12	医薬品等原薬	メチルプレドニゾロンコハク酸エステルナトリウム	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	2375-03-3	CC1CC2C3CCC(C3(CC(C2C4(C1=CC(=O)C=C4)C)O)C)(C(...	<img data-content="rdkit/molecule" src="data:i...
52	220MF10128	2008/05/15	2019/03/05	医薬品等原薬	メチルプレドニゾロンコハク酸エステル	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	False	False	False
...	...	...	...	...	...	...	...	...	...	...	...	...
2633	219MF10105	2007/03/23	2012/03/30	医薬品等原薬	イトラコナゾール	ｍｙｌａｎ　ｌａｂｏｒａｔｏｒｉｅｓ　ｌｉｍｉｔｅｄ	ｐｌｏｔ　ｎｏ．５６４／ａ／２２，　ｒｏａｄ　ｎｏ．９２，　ｊｕｂｉｌｅｅ　ｈｉｌｌｓ，　ｈ...	NaN	NaN	False	False	False
2821	222MF10253	2010/10/28	2010/10/28	医薬品等原薬	オロパタジン塩酸塩	ｋｏｌｏｎ　ｌｉｆｅ　ｓｃｉｅｎｃｅ，ｉｎｃ．	１１０，ｍａｇｏｋｄｏｎｇ－ｒｏ，ｇａｎｇｓｅｏ－ｇｕ，ｓｅｏｕｌ，ｋｏｒｅａ	NaN	NaN	False	False	False
3478	217MF10151	2005/05/31	2008/08/18	医薬品等原薬	マレイン酸セチプチリン「サンヨー」	株式会社　三洋化学研究所	大阪府堺市美原区多治井１４８－１	NaN	NaN	False	False	False
3793	219MF10159	2007/05/15	2007/05/15	医薬品等原薬	チクロピジン塩酸塩「テバ」	ｓｉｃｏｒ　ｓｏｃｉｅｔａ（｀）　ｉｔａｌｉａｎａ　ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　ｓ．ｒ．ｌ．	ｐｉａｚｚａｌｅ　ｌｕｉｇｉ　ｃａｄｏｒｎａ，　４　－　２０１２３　ｍｉｌａｎｏ，　ｉｔａｌｙ	NaN	NaN	False	False	False
3809	218MF10270	2006/02/17	2007/04/27	医薬品等原薬	日本薬局方　メキタジン	住友化学株式会社	東京都中央区日本橋二丁目７番１号	NaN	NaN	False	False	False
1129 rows × 12 columns
"""

最後に、こちらを参考にして合成難易度 (sa)を計算しました。MolObjectを受け取ってsaを返す関数rtn_saを作成しました。

def rtn_sa(mol):
    if mol==False:
        return False
    sa = sascorer.calculateScore(mol)
    return sa

df_test.loc[:,'sa'] = df_test['ROMol'].map(rtn_sa)
df_test

"""
	MF登録番号	初回登録年月日	最新登録年月日	登録区分	登録品目名	登録者氏名	登録者住所	原薬等国内管理人の名称	原薬等国内管理人の住所	cas	smiles	ROMol	sa
26	219MF10293	2007/09/10	2019/03/29	医薬品等原薬	スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	69388-84-7	CC1(C(N2C(S1(=O)=O)CC2=O)C(=O)[O-])C.[Na+]	<img data-content="rdkit/molecule" src="data:i...	4.448697
28	218MF10710	2006/08/14	2019/03/29	医薬品等原薬	日本薬局方　スルバクタムナトリウム（製造専用）	ｋｙｏｎｇｂｏ　ｐｈａｒｍａｃｅｕｔｉｃａｌ　ｃｏ．，ｌｔｄ．	１７４，ｓｉｒｏｋ－ｒｏ，ａｓａｎ－ｓｉ，ｃｈｕｎｇｃｈｅｏｎｇｎａｍ－ｄｏ，ｋｏｒｅａ	NaN	NaN	False	False	False	False
36	218MF10971	2006/12/01	2019/03/19	医薬品等原薬	ナロキソン塩酸塩水和物	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	51481-60-8	C=CCN1CCC23C4C(=O)CCC2(C1CC5=C3C(=C(C=C5)O)O4)...	<img data-content="rdkit/molecule" src="data:i...	5.410687
42	230MF10133	2018/09/14	2019/03/12	医薬品等原薬	メチルプレドニゾロンコハク酸エステルナトリウム	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	2375-03-3	CC1CC2C3CCC(C3(CC(C2C4(C1=CC(=O)C=C4)C)O)C)(C(...	<img data-content="rdkit/molecule" src="data:i...	5.00208
52	220MF10128	2008/05/15	2019/03/05	医薬品等原薬	メチルプレドニゾロンコハク酸エステル	ｓａｎｏｆｉ　ｃｈｉｍｉｅ	８２　ａｖｅｎｕｅ　ｒａｓｐａｉｌ　９４２５０　ｇｅｎｔｉｌｌｙ　ｆｒａｎｃｅ	NaN	NaN	False	False	False	False
...	...	...	...	...	...	...	...	...	...	...	...	...	...
2633	219MF10105	2007/03/23	2012/03/30	医薬品等原薬	イトラコナゾール	ｍｙｌａｎ　ｌａｂｏｒａｔｏｒｉｅｓ　ｌｉｍｉｔｅｄ	ｐｌｏｔ　ｎｏ．５６４／ａ／２２，　ｒｏａｄ　ｎｏ．９２，　ｊｕｂｉｌｅｅ　ｈｉｌｌｓ，　ｈ...	NaN	NaN	False	False	False	False
2821	222MF10253	2010/10/28	2010/10/28	医薬品等原薬	オロパタジン塩酸塩	ｋｏｌｏｎ　ｌｉｆｅ　ｓｃｉｅｎｃｅ，ｉｎｃ．	１１０，ｍａｇｏｋｄｏｎｇ－ｒｏ，ｇａｎｇｓｅｏ－ｇｕ，ｓｅｏｕｌ，ｋｏｒｅａ	NaN	NaN	False	False	False	False
3478	217MF10151	2005/05/31	2008/08/18	医薬品等原薬	マレイン酸セチプチリン「サンヨー」	株式会社　三洋化学研究所	大阪府堺市美原区多治井１４８－１	NaN	NaN	False	False	False	False
3793	219MF10159	2007/05/15	2007/05/15	医薬品等原薬	チクロピジン塩酸塩「テバ」	ｓｉｃｏｒ　ｓｏｃｉｅｔａ（｀）　ｉｔａｌｉａｎａ　ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ　ｓ．ｒ．ｌ．	ｐｉａｚｚａｌｅ　ｌｕｉｇｉ　ｃａｄｏｒｎａ，　４　－　２０１２３　ｍｉｌａｎｏ，　ｉｔａｌｙ	NaN	NaN	False	False	False	False
3809	218MF10270	2006/02/17	2007/04/27	医薬品等原薬	日本薬局方　メキタジン	住友化学株式会社	東京都中央区日本橋二丁目７番１号	NaN	NaN	False	False	False	False
1129 rows × 13 columns
"""

データ分析するために前処理を行います。saが計算できなかったものは除外しました。Top20で抽出しましたが、処理の過程で一社欠損してしまったようでした。

df_ana  = df_test[df_test['sa']!=False]
print(df_ana.shape)
print(len(df_ana['登録者氏名'].unique()))
"""
(408, 13)
19
"""

saについて統計量を計算しました。平均値についてランキングして整理した結果が以下の通りです。社名に\u3000が入り込んでいるのは気にしないでください。

df_rank = df_ana.groupby(['登録者氏名'])['sa'].agg(['mean', 'std', 'count', 'max', 'min']).sort_values(by="mean", ascending=False)
df_rank

	mean	std	count	max	min
登録者氏名
ｓａｎｏｆｉ\u3000ｃｈｉｍｉｅ	4.791725	1.956721	22	8.378029	2.243329
アルプス薬品工業株式会社	4.613620	0.906464	15	5.929040	2.869499
ｈａｎｍｉ\u3000ｆｉｎｅ\u3000ｃｈｅｍｉｃａｌ\u3000ｃｏ．，ｌｔｄ．	3.917504	1.021601	14	6.116182	2.505240
ｋｙｏｎｇｂｏ\u3000ｐｈａｒｍａｃｅｕｔｉｃａｌ\u3000ｃｏ．，ｌｔｄ．	3.655462	1.243465	39	8.256955	1.871614
ｋｏｌｏｎ\u3000ｌｉｆｅ\u3000ｓｃｉｅｎｃｅ，ｉｎｃ．	3.620023	0.987406	16	5.361381	2.492740
ｓｉｃｏｒ\u3000ｓｏｃｉｅｔａ（｀）\u3000ｉｔａｌｉａｎａ\u3000ｃｏｒｔｉｃｏｓｔｅｒｏｉｄｉ\u3000ｓ．ｒ．ｌ．	3.347663	1.118844	4	4.828636	2.122733
ｈｅｔｅｒｏ\u3000ｌａｂｓ\u3000ｌｉｍｉｔｅｄ	3.278368	0.845611	22	5.326849	2.337358
ｍｙｌａｎ\u3000ｌａｂｏｒａｔｏｒｉｅｓ\u3000ｌｉｍｉｔｅｄ	3.259816	0.926031	32	6.116182	2.251534
相模化成工業株式会社	3.246741	1.018134	13	5.086630	2.116077
ａｕｒｏｂｉｎｄｏ\u3000ｐｈａｒｍａ\u3000ｌｉｍｉｔｅｄ	3.207911	0.735930	19	4.445804	2.343852
株式会社\u3000大阪合成有機化学研究所	3.205708	1.046548	12	5.629029	2.257189
ｏｌｏｎ\u3000ｓ．ｐ．ａ．	3.119780	0.753195	24	5.625602	2.275611
ダイト株式会社	3.053274	0.791519	30	5.361381	2.097384
株式会社\u3000三洋化学研究所	3.007404	0.630739	41	5.361381	2.116077
株式会社パーマケム・アジア	2.914877	0.701407	31	4.794299	1.469664
ｃａｍｂｒｅｘ\u3000ｐｒｏｆａｒｍａｃｏ\u3000ｍｉｌａｎｏ\u3000ｓ．ｒ．ｌ．	2.850941	1.009620	23	7.173955	2.113785
住友化学株式会社	2.741267	0.498912	19	3.580353	1.811313
大和薬品工業株式会社	2.733274	0.440762	11	3.342106	1.819255
八代製薬株式会社	2.590104	0.610828	21	3.705913	1.407299

Top20中のsaの分布を以下に示します。sa=3を境にグッと分布が少なくなっています。

plt.grid()
plt.hist(df_ana['sa'], bins=20)
plt.xlim(0,10)
plt.xlabel('sa')
plt.ylabel('count')
plt.xticks(np.arange(0,10,1))
plt.title('Distribution of total sa')

最後に、企業別のsaの分布を以下に示します。タイトルが長いと別のグラフに被ってしまうので、ある程度の長さで切っています。

fig, axes = plt.subplots(4,5, figsize=(15,15))

for com, ax in zip(df_ana['登録者氏名'].unique(), axes.ravel()):
    df_tmp = df_ana[df_ana['登録者氏名']==com]
    ax.hist(df_tmp['sa'].astype(float).round(2), bins=10)
    ax.set_xlim(0,10)
    ax.set_ylim(0,12)
    ax.set_xlabel('sa')
    ax.set_ylabel('count')
    ax.set_xticks(np.arange(0,10,1))
    ax.set_yticks(np.arange(0,12,1))
    if len(com)>=15:
        com = com[:15]
    ax.set_title(com)
    ax.grid()
fig.tight_layout()

考察

今回の調査で得た考察を以下に示します。

• 全てのデータからTop20の企業を抽出すると、1129件であった。そこから、登録品目名+casというキーワードでgoogle検索してcas番号を正規表現で取り出した。検索画面からcas番号を見つけ出すため、見逃しが多少あった。cas番号からsmilesの取得をpubchemで行った。この検索もヒットしない場合もあったため、最終的には408件のデータを抽出した。
• Synthetic Accessibility score (sa：高いほど合成難易度が高い)について、企業ごとに基本統計量 (平均、標準偏差、個数、最大値、最小値)を求めた。ｓｉｃｏｒ　ｓｏｃｉｅｔａはn=4と最もエントリー数が少なかったが、その他のエントリー数は10個以上であったので、ある程度意味のある比較が可能であると判断した。
• 全体のsaの分布を確認すると、およそ3が最頻値であり、それ以上の値は極端に少なくなっていたことから、sa=3が一つの指標となる。
• saの平均値でランキングしたところ、1位はsanofiであり、最大値も1位であった。エントリー数も22個であり、sa=8の化合物が4件もあることから、間違いなく技術力が高い企業であることが読み取れる。
• MF申請数ランキング実質1位であった三洋化学研究所は、saランキングでは14位となっており、あまりチャレンジせずに手堅く実績を残している様子であった。
• MF申請数ランキング実質3位であったアルプス薬品工業は、saランキングでも2位と海外勢を抑えて上位に位置しており、質と量ともにレベルの高さがうかがえる。標準偏差も0.9と常に挑戦的な研究開発に取り組んでいると考えられる。
• MF申請数ランキング実質5位であったダイトは、saランキングでは13位となっており、あまりチャレンジせずに手堅く実績を残している様子であった。
• 上位の半数以上が海外勢であり、上記のバブルチャートでも述べた通り、研究開発は世界レベルのものが求められている。

所感

今回のデータ分析の結果から、国内の競合他社は三洋化学研究所、アルプス薬品工業、ダイトであることが分かりました。国外からの参入もあり、勝ち残っていくには世界レベルの研究開発が必要であることも分かりました。
コード技術の観点からでは、表記されている化合物から合成難易度を算出して、企業の研究開発レベルを定量化する方法は我ながらよく思いついたなと感じました。企業情報をさらに付け加えたいと思ったのですが、企業データベースは基本有料だったので、止む無く断念しました。

今後の予定

医薬品業界を知る情報源はまだまだありますので、今後も情報発信をしていきたいと思います。
最後まで読んでいただきありがとうございました。

参考文献

【特許分析】Pythonでお金を掛けずにパテントマップを作ってみた
RDKitで合成難易度を評価して化合物をスクリーニング (化学の新しいカタチ)