ニューラルネットワークによる多種類の有機化合物の発ガン性の予測

田辺 和俊, 大森 紀人, 小野 修一郎, 鈴木 孝弘, 松本 高利, 長嶋 雲兵, 上坂 博亨

1 緒 言

2 方法

2. 1 ニューラルネットワークの構造

2. 2 発ガン性データ

Table 1. Numbers of compounds

Data Set	Carcinogen	Non-carcinogen	Total
Model construction	107	180	287
Model evaluation	45	122	167
Total	152	302	454

2. 3 記述子データ

Table 2. Descriptors used in this study

Descriptor Set	Descriptors	Ra
tReymers	SLOGP	-0.01
	SMR	-0.03
	ROTBOND	-0.06
	FLEX	0.02
	VOLUME	-0.09
	SURF_A	-0.09
	TPSA	-0.08
	HBD	-0.06
	HBA	-0.04
	LUMO	0.03
	HOMO	-0.01
	DIPOLE	-0.15
	LOGD_2	-0.05
	LOGD_7	-0.02
	LOGD_10	0.00
Helma	logP	0.02
	dipole	-0.08
	electronic energy	0.06
	electronegativity	-0.04
	heat of formation	0.17
	homo	-0.06
	homo-lumo	-0.04
	hybrid dipole	-0.01
	ionization potential	0.06
	lumo	0.00
	largest interatomic distance	-0.07
	molecular weight	0.02
	point-chg. dipole	-0.09
	total energy	0.05
	/POLARIZA	-0.04
	/DELTAHF	0.19
	/STABIL	0.03
	/STRAIN	-0.11
	total_access	-0.07
	non_polar_access	-0.04
	polar_access	-0.04
	perc_nonpolar	0.04

^aR : Correlation coefficient with carcinogenicity

QSARによる予測モデルを構築する際、最も重要な因子は入力記述子の選択である。Table 2に示すように、今回用いたPTCの記述子の中にはhomo、lumo、log PのようにtReymersとHelmaで重複していたり、相関が高いと思われる記述子が幾つかある。しかし、それらの記述子はTable 3の相関係数に示すように、tReymersとHelmaでの値が異なっている。

Table 3. Correlation coefficients between similar descriptors

Descriptor A	Descriptor B	Ra
tReymers SLOGP	Helma logP	0.94
tReymers SMR	Helma /POLARIZA	1.00
tReymers SURF_A	Helma total_access	0.94
tReymers LUMO	Helma lumo	0.14
tReymers HOMO	Helma homo	0.76
tReymers DIPOLE	Helma dipole	0.76
Helma heat of formation	Helma /DELTAHF	0.90
Helma homo	Helma ionization potential	-1.00

^aR: Correlation coefficient between descriptors A and B

これはそれらの記述子を計算した時の化合物の立体配座 (構造最適化の計算方法の違いによる) などが異なるためと思われるが、我々としてどちらの値を用いるべきかの選択にあたっての基準が見当たらない。そこで、本研究ではこのような重複している記述子から独立かつ有効な情報量を抽出するために、主成分分析を行い、その結果の累積寄与率から判断した主成分の値をニューラルネットワークの入力層に入力した。
ただし、主成分分析の結果をニューラルネットワークの入力に使うと、記述子の感度分析が行えなくなるなどの問題が生じる。しかし、今回用いた37個の記述子の間にはかなり線形の相関が高い記述子の組がいくつかあり、また学習用のデータについて自由度（すなわち、データの数とパラメータ数の比）を高く保つため、入力変数の数を絞り込む必要があった。そこで、本研究では主成分分析を行うことにしたが、その有効性を確認するために、主成分分析を行わなかった結果と比較した。
一方、出力層の教師データには発ガン性化合物に0.9、非発ガン性化合物に0.1の値を入力した。このように、教師データを0.0 ～ 1.0でなく、0.1 ～ 0.9に規格化したのは、ニューラルネットワークのユニット (ニューロン) において数値変換に使われているsigmoid関数の特性を考慮したためである。また、テスト化合物の予測的中率の判定においては、出力層の値が0.5以上なら発ガン性あり、0.5以下なら発ガン性なしと判定した。

2. 4 ニューラルネットワークの学習

3 結果と考察

3. 1 入力記述子の選定

Table 4. Result of principal component analysis

No.	Eigenvalue	Proportion (%)	Cumulative proportion (%)
1	0.2284	36.0084	36.0084
2	0.1356	21.3767	57.3851
3	0.0819	12.9115	70.2966
4	0.0469	7.3974	77.6941
5	0.0300	4.7265	82.4206
6	0.0227	3.5843	86.0049
7	0.0135	2.1282	88.1331
8	0.0115	1.8156	89.9487
9	0.0103	1.6201	91.5688
10	0.0096	1.5073	93.0761
11	0.0073	1.1503	94.2264
12	0.0060	0.9441	95.1705
13	0.0044	0.6990	95.8696
14	0.0037	0.5794	96.4490
15	0.0032	0.5081	96.9571

3. 2 学習回数の決定

Figure 1. Variations of average errors of training (bold line), validation (thick line) and test (thin line) sets in training cycles

この図から、error-back-propagation法による学習を繰り返すと、学習用化合物の誤差は限りなく減少していくのに対し、検証用化合物の誤差は最初は減少していくが、ある回数において極小を示し、その後は増加に転じていることがわかる。この検証用化合物の誤差が増加した状況はニューラルネットワークの過学習を示しており、この誤差が極小を示す回数が最適の学習回数であると考えられる。本研究ではこの方法により、ニューラルネットワークの過学習を解決し、最適の学習回数を決定した。

3. 3 中間層ユニット数の決定

Figure 2. Dependence of average error of validation set on numbers of intermediate layer units

この問題に対して栗田 [18] はAIC (Akaike Information Criteria) 理論 [19]を適用して中間層ユニット数の最適値を決定する方法を提案している。本研究の3層構造のニューラルネットワークの場合には、AICは次式のように定義される：
AIC = Nk log(σ²) + 2K
ここで、Nは学習データの数、kは出力層のユニット数、σ²は出力値と教師データとの誤差の分散、Kはニューラルネットワークに含まれるパラメータの数である。今の場合に上記の誤差の結果を用いてAICを計算し、中間層ユニット数に対してプロットするとFigure 3のようになった。この結果からAICが最小になる中間層ユニット数の最適値を2に決定した。

Figure 3. Dependence of AIC on numbers of intermediate layer units

3. 4 重み初期値の決定

Table 5. Results of our model^a

	RMSEb			CCRc of Training Set			CCRc of Validation Set			CCRc of Test Set
Random Seed	Cycled	Tre	Vae	Tee	+f	-f	Overall	+f	-f	Overall	+f	-f	Overall
9199	146	0.1058	0.1436	0.1306	0.113	1.000	0.671	0.056	0.978	0.632	0.044	0.992	0.737
9116	199	0.1058	0.1439	0.1296	0.151	0.967	0.664	0.056	0.956	0.618	0.000	0.992	0.725
16931	250	0.1072	0.1439	0.1336	0.132	1.000	0.678	0.111	0.967	0.646	0.044	0.992	0.737
22801	179	0.1055	0.1440	0.1305	0.094	1.000	0.664	0.093	0.989	0.653	0.044	0.992	0.737
10111	2604	0.1058	0.1440	0.1747	0.491	0.856	0.720	0.389	0.789	0.639	0.489	0.623	0.587
22157	277	0.1055	0.1441	0.1339	0.075	1.000	0.657	0.037	0.967	0.618	0.067	0.992	0.743
18382	380	0.1082	0.1441	0.1303	0.094	1.000	0.664	0.037	0.978	0.625	0.000	0.992	0.725
2417	441	0.1053	0.1441	0.1292	0.132	0.956	0.650	0.093	0.956	0.632	0.000	0.992	0.725
4418	154	0.1068	0.1441	0.1315	0.094	0.989	0.657	0.037	0.978	0.625	0.000	1.000	0.731
29762	258	0.1138	0.1441	0.1301	0.094	1.000	0.664	0.130	0.967	0.653	0.000	0.975	0.713

^aUnit number of intermediate layer=2^bRMSE : Root mean square error between output and teaching data^cCCR : Correct classification rate^dOptimum training cycle where RMSE of the validation set was minimum.^eTr: training set, Va: validation set, Te: test set^f+ : carcinogenic active, - : carcinogenic inactive

3. 5 結果の比較

Table 6. Correct classification rates of various models

No	Model name	Carcinogen	Non-carcinogen	Overall
1	This work	0.044	0.992	0.737
2	ANI1	0.100	0.870	0.664
3	ANI2	0.077	0.941	0.710
4	ANI3	0.159	0.860	0.672
5	ANU1	0.611	0.402	0.458
6	ANU2	0.350	0.653	0.572
7	ANU3	0.572	0.300	0.373
8	BAUS	0.312	0.773	0.649
9	BAUN	0.251	0.771	0.631
10	BAUP	0.352	0.712	0.616
11	GONS	0.251	0.850	0.689
12	JAN0	0.363	0.671	0.589
13	JAN4	0.240	0.732	0.601
14	KWA1	0.169	0.909	0.711
15	LEU1	0.900	0.071	0.293
16	LEU2	0.098	0.832	0.635
17	LISI	0.271	0.740	0.615
18	LUC1	0.503	0.501	0.502
19	LUC2	0.210	0.801	0.643
20	LUC3	0.210	0.771	0.621
21	PLE1	0.253	0.633	0.531
22	PLE2	0.141	0.781	0.609
23	PLE3	0.391	0.463	0.443
24	SMU1	0.350	0.653	0.572
25	SMU2	0.210	0.740	0.598
26	SMU3	0.523	0.483	0.494
27	SUT	0.523	0.531	0.529
28	VINI	0.411	0.663	0.596
29	WAI1	0.652	0.412	0.476
30	WAI2	0.332	0.694	0.597
31	WAI3	0.442	0.552	0.522

ただし、上記のように、本研究の的中率は、テスト用化合物について最適化したモデルで得られたのではないという点では、コンテストの参加者と同じではあるが、参加者は正解を知らない条件で解析しているので、比較は公平ではない。参加者のモデルの結果のうちで公表されていない結果が本研究の的中率を超える可能性があることは否定できない。
本研究の的中率が高くなった理由としては、記述子の選択とニューラルネットワークという非線形解析モデルの利用の2点が考えられる。参加者の予測モデルの大半では、前記のように非常に多数の記述子と比較的単純な線形解析モデルが用いられているため、構造情報のみから多種類の化合物の発ガン性を予測するという問題では、このような機械的な解析モデルでは良好な成績を得にくいと考えられる。
すなわち、上記のようにコンテストの主催者から公開された7000種類以上にのぼる記述子の大半は、化合物中の種々の原子、結合、部分構造の有無、個数などを示す記述子であり、ほとんど0の値をもつsparseな記述子である。このような記述子を用いて比較的単純な線形モデルで機械的に用いて解析しても良好な成績を得ることは困難であろう。
それに対し、本研究では記述子の内容を吟味し、このようなsparseな記述子は排除し、Table 2に示す記述子のみを用いて解析を行い、参加者より良好な結果を得た。このような記述子選択の重要性については、今回のコンテスト以前に開発されていた既存のモデル、たとえばCASE [20 - 24]、TOPKAT [25 - 28]、REPAD [29]、COMPACT [30]、FALS [31 - 33]などの性能が十分でない原因として、原子、結合、部分構造の数などのsparseな記述子が用いられていることからも支持される。

3. 6 発ガン性データの問題

Figure 4. ROC points and convex hull of various models (Numbers in the figure correspond to those in Table 6)

しかし、我々のモデルおよび上記の2モデルでは、テスト化合物のうちで非発ガン性化合物の的中率はいずれも好成績であるのに対し、発ガン性化合物の的中率はきわめて低い。このような観点から、Table 6およびFigure 4の結果を分析すると、我々とコンテスト参加者の成績は大まかに3つのグループに分けられる。
第1グループは、我々及びANI、KWAであり、上記のように発ガン予測率はきわめて低いが、非発ガン予測率が高く、両者を合わせた予測率も高いグループである。第2グループは、ANU、LEU、WAIであり、第1グループとは逆に、発ガン予測率は高いが、非発ガン予測率がきわめて低く、両者を合わせた予測率もかなり低いグループである。第3グループは、残りのモデルで、発ガン予測率も非発ガン予測率も両者を合わせた予測率も中程度のグループである。このような分類は、各モデルの記述子の選択あるいは相関解析法の違いを反映していると考えられる。
我々のモデルの場合、学習用、検証用、テスト用各化合物に対するニューラルネットワークの出力値の分布をFigure 5に示す。この図のテスト用化合物の出力値の分布から、非発ガン性化合物は好成績で予測できるのに対し、発ガン性化合物がほとんど予測できていないことが分かる。また、PTCコンテストのその他の大半のモデルでも発ガン性化合物より非発ガン性化合物の的中率が高い。

Figure 5. Histogram of output values of neural network for carcinogen ( ) and non-carcinogen ( ) of training, validation and test sets

この結果は実用的にはきわめて問題である。すなわち、冒頭に記したように、QSARモデルによる有害性予測を動物試験の前段階のスクリーニングと位置付けるならば、発ガン性化合物を非発ガンと判定して見逃すことは好ましくないからである。すなわち、非発ガン性化合物よりも発ガン性化合物の方をできるだけ高い精度で予測することが実用上からは望ましいからである。
このように、PTCコンテストのデータを用いた我々のモデルを含む大半のモデルが、非発ガン性は高成績で、しかし発ガン性は低成績でしか予測できなかった原因としては、Table 1に示すように、モデル構築用化合物の数が発ガン性107 種類、非発ガン性180種類と不均衡なために、モデルの学習が非発ガン性化合物に偏って行われたためと考えられる。
また、Table 6に見られるように、PTCコンテストの応募者の中には発ガン性化合物の的中率が非発ガン性化合物より高くなったモデルも幾つかある。しかし、それらのモデルによる発ガン性、非発ガン性、両者を合計した的中率はいずれも50%前後であり、QSARによる発ガン性予測モデルとしてはきわめて不十分な成績である。この原因としては、NTPのモデル構築用化合物群が構造的に広範囲なバルクの化学物質を含み、かつ比較的小さな分子が多いのに対し、FDAのテスト化合物群は医薬品に偏り、化学構造の面でかなり異なっている点が指摘されている [9, 10]。

3. 7 記述子の問題

3. 8 今後の課題

4 結論

本研究は科学研究費補助金基盤研究 (A)(1) 14209022により行われた。また、情報量基準について助言いただいた産業技術総合研究所脳神経情報研究部門の栗田多喜夫博士に感謝します。

参考文献