div.ProseMirror

DeepChem 라이브러리를 활용한Solubility Modeling

version check

!python --version
0.9s
Python
→
import tensorflow as tf
tf.__version__
3.4s
Python
→
'1.14.0'
import deepchem as dc
dc.__version__
2.6s
Python
→
'2.3.0'

Modeling Solubility

분자 용해도를 계산하여 예측하는 것은 약물 발견에 유용합니다. 약물과 유사한 화합물의 용해도를 예측하는 모델을 구축하는데 deepchem 라이브러리를 사용합니다. 모델에 fitting하는 과정에는 4단계가 있습니다.

  1. 수성 용해도 측정과 함께 일련의 화합물로 구성된 화학 데이터세트 로드

  2. 각 화합물을 통계 학습 방법이 이해할 수 있는 특징 벡터(inline_formula not implemented)로 변환

  3. 특징 벡터를 수용성 용해도의 추정치에 매핑하는 간단한 모델 학습

  4. 결과 시각화

from rdkit import Chem
from rdkit.Chem import Draw
from rdkit.Chem.Draw import IPythonConsole
from rdkit.Chem import Descriptors
from rdkit.Chem import AllChem
from rdkit import DataStructs
import numpy as np
0.2s
Python
→

Load to data

우리는 추정된 수용성 용해도 측정 데이터 세트[1]를 Deepchem에 로드해야합니다. 데이터는 csv 형식이고 SMILES 문자열, 예측된 수성 용해도 및 많은 분자특성을 포함합니다.

이 필드의 대부분은 우리의 목적에 유용하지 않습니다. 우리가 필요로하는 두가지 속성은 SMILEs 와 Measured log-solubility in mols/liter입니다. 데이터를 Deepchem에 로드하기 전에 데이터를 python에 로드하고 간단한 예비 분석을 수행하여 데이터에 대한 이해를 해봅시다.

!wget https://raw.githubusercontent.com/deepchem/deepchem/master/datasets/delaney-processed.csv
#!wget -O delaney-processed.csv https://raw.githubusercontent.com/deepchem/deepchem/master/datasets/delaney-processed.csv
#!wget -P /want/directory https://raw.githubusercontent.com/deepchem/deepchem/master/datasets/delaney-processed.csv
1.1s
Python
→
#데이터가 잘 다운로드되었는지 확인
!ls
0.9s
Python
→
from deepchem.utils.save import load_from_disk
dataset_file= "delaney-processed.csv"
dataset = load_from_disk(dataset_file)
print("Columns of dataset: %s" % str(dataset.columns.values))
print("Number of examples in dataset: %s" % str(dataset.shape[0]))
0.4s
Python
→
print(type(dataset))
dataset
0.5s
Python
→

Data EDA

데이터 세트에서 화합물을 시각적으로 이해하려면 RDkit을 사용하여 화합물을 그려 봅시다.

우선 이미지로 출력하기 위한 함수를 정의가 아래 코드에 있습니다.

import tempfile
from rdkit import Chem
from rdkit.Chem import Draw
from itertools import islice
from IPython.display import Image, display
# 이미지로 출력하기 위한 함수
def display_images(filenames):
    """Helper to pretty-print images."""
    for file in filenames:
      display(Image(file))
#SMILES데이터를 Mol Object로 변환 후 mol Object list를 입력으로 주면 해당 리스트의 분자를 그림으로 출력
def mols_to_pngs(mols, basename="test"):
    """Helper to write RDKit mols to png files."""
    filenames = []
    for i, mol in enumerate(mols):
        filename = "%s%d.png" % (basename, i)
        Draw.MolToFile(mol, filename)
        filenames.append(filename)
    return filenames
0.0s
Python
→

입력 받은 데이터의 분자를 시각화해봅시다.

num_to_display = 5
molecules = []
for _, data in islice(dataset.iterrows(), num_to_display):
    molecules.append(Chem.MolFromSmiles(data["smiles"]))
display_images(mols_to_pngs(molecules))
0.4s
Python
→

입력 받은 데이터의 타겟변수 logS 분포를 살펴봅시다. 내가 사용하고자하는 함수의 기능이나 인자를 살펴보기위해서는 인터넷에 패키지명과 함수명으로 검색을 수행하거나, 코드 작업 중인 창에서 함수명 뒤에 커서를 올려놓고 기다리면 함수 doc을 볼 수 있습니다.

import matplotlib
import numpy as np
import matplotlib.pyplot as plt
solubilities = np.array(dataset["measured log solubility in mols per litre"])
#히스토그램을 그리기 위한 함수 사용 단일 변수를 입력으로 주면 분포 히스토그램을 출력해주는 기능
#인터넷에 패키지 이름과 함수명으로 검색하면 아래와 같은 doc을 찾을 수 있음
# https://matplotlib.org/3.3.1/api/_as_gen/matplotlib.pyplot.hist.html
n, bins, patches = plt.hist(solubilities, 50, facecolor='green', alpha=0.75)
#출력할 그래프의 x축 y축 이름을 설정하는 부분
plt.xlabel('Measured log-solubility in mols/liter')
plt.ylabel('Number of compounds')
#그래프의 제목
plt.title('Histogram of solubilities')
#격자 무늬 출력
plt.grid(True)
#그래프 출력 명령어
plt.show()
0.8s
Python
→
import seaborn as sns
from scipy.stats import norm
solubilities = np.array(dataset["measured log solubility in mols per litre"])
#히스토그램을 그리기 위한 함수 사용 단일 변수를 입력으로 주면 분포 히스토그램을 출력해주는 기능
sns.distplot(solubilities, color="blue", kde=True)
plt.xlabel('Measured log-solubility in mols/liter')
plt.ylabel('Number of compounds')
#그래프의 제목
plt.title('Histogram of solubilities')
#격자 무늬 출력
plt.grid(True)
#그래프 출력 명령어
plt.show()
0.6s
Python
→

데이터 변환(Transform)

사전 분석이 완료되면, deepChem라이브러리를 사용하여 분자 용해도의 예측 모델을 구성합니다. 이러한 분자를 생성하는 첫 번째 단게는 각 화합물을 통계적 학습 기술로 이해할 수 있는 벡터 형식으로 변환하는 것입니다. 이러한 과정을 일반적으로 featurization(특성화)라고 합니다.

대부분의 얻을 수 있는 데이터에는 분자구조 정보가 SMILES형식으로 되어있습니다. 이에 대한 특성 벡터를 생성하는 방법으로 ECPF를 사용합니다.[3] deepChem은 특성화를 위한 API를 제공하고 이는 dc.feat.하위에 대부분 존재합니다. 아래 코드의 순서는 다음과 같습니다.

  1. 특성화를 시작하기 위해 먼저 특성화 객체를 생성합니다. deepChem은 CircularFingerprint 클래스(ECFP 특성화를 수행하는 Featurizer의 하위 클래스)를 제공합니다.

  2. 데이터를 불러오고 불러온 데이터에 특성화 작업을 수행합니다.

import deepchem as dc
featurizer = dc.feat.CircularFingerprint(size=1024)
0.0s
Python
→

이제 실제 특성화를 수행 해봅시다. deepChem은 이러한 목적으로 CSVLoader 클래스를 지원합니다. 이 클래스의 featurize()메서드는 디스크에서 데이터를 로드하고 제공된 Featurizer 인스턴스를 사용하여 제공된 데이터를 특징 벡터로 변환합니다.

이러한 데이터 세트에 대해 기계학습을 수행하려면 샘플을 기계 학습에 적합한 데이터세트로 변환해야합니다.(즉, 데이터 행렬 inline_formula not implemented, 여기서 inline_formula not implemented은 샘플 수를 inline_formula not implemented는 특징 벡터의 차원을 레이블 벡터는 다음과 같이 표현합니다. inline_formula not implemented)

#첫번째 인자는 해당 SMILES의 Label로 사용할 열의 이름
#두번째 인자는 변환할 SMILES 데이터가 존재하는 열의 이름
#세번째 인자는 앞서 생성한 Featurizer 객체
loader = dc.data.CSVLoader(
      tasks=["measured log solubility in mols per litre"], smiles_field="smiles",
      featurizer=featurizer)
dataset = loader.featurize(dataset_file)
6.4s
Python
→
print(dataset.X[0], dataset.y[0])
print(len(dataset.X[0]), len(dataset.y[0]))
print(type(dataset))
0.3s
Python
→

데이터 분할(Split Data)

모델을 구성할 때는 제공된 데이터를 Training/Test subset으로 분리해야합니다.

Training set은 모델을 학습하는데 사용되고 Test set은 학습된 모델을 평가하는데 사용됩니다.

실제로 모델 학습에는 그림 하단부와 같이 Training, Validation, Test split을 수행하는 것이 대체로 유용합니다.

모델의 학습과 학습과 동시에 모델 성능 평가의 에는 Training set과 Validation set을 사용합니다. Testset은 주어진 데이터에서 Split하고, 어떤 데이터인지 우리는 알고 있지만, 학습하는 모델이 모르게 해야합니다.

훈련, 검증, 테스트 분할을 수행하는 적절한 방법을 선택하는 것은 어려울 수 있습니다.

기계학습 모델의 입력 데이터 분할 시 데이터를 Train/Validation/Test로 Random split하는 것이 일반적 입니다. 하지만, 무작위 분할은 화학 정보학의 목적에 알맞지 않습니다.

예측 모델이 유용하기 위해서는 훈련 데이터의 분자 세트를 넘어서 화학 공간의 일부에 예측력이 있어야 합니다.

결론적으로, 데이터 분할 과정에서 Bemis-Murcko Scaffolds[5]를 사용하여 이 분할을 수행합니다. 기본 분자 구조(scaffolds)를 공유하는 모든 화합물은 Train/Validation/Test split에서 동일한 분할 위치에 배치됩니다. 꼭 Scaffolds Split이 유용한 것은 아니지만 이번에는 Scaffolds Split을 사용합니다.

또한 Moleculenet의 benchmark 지침을 따릅니다.

splitter = dc.splits.ScaffoldSplitter()
train_dataset, valid_dataset, test_dataset = splitter.train_valid_test_split(
    dataset)
0.8s
Python
→

정말로 train과 test split의 분자 구성이 다른지 살펴봅시다.

train_dataset.ids[:5]
0.0s
Python
→
array(['CC(C)=CCCC(C)=CC(=O)', 'c2ccc1scnc1c2 ', 'Clc1cc(Cl)c(c(Cl)c1)c2c(Cl)cccc2Cl', 'CC12CCC3C(CCc4cc(O)ccc34)C2CCC1O', 'Clc1ccc2ccccc2c1'], dtype=object)
train_mols = [Chem.MolFromSmiles(compound)
              for compound in train_dataset.ids]
display_images(mols_to_pngs(train_mols[:5], basename="train"))
1.3s
Python
→
test_dataset.ids[:5]
0.2s
Python
→
array(['OCC3OC(OCC2OC(OC(C#N)c1ccccc1)C(O)C(O)C2O)C(O)C(O)C3O ', 'Cc1occc1C(=O)Nc2ccccc2', 'c1ccc2c(c1)ccc3c2ccc4c5ccccc5ccc43', 'c1ccsc1', 'ClC4=C(Cl)C5(Cl)C3C1CC(C2OC12)C3C4(Cl)C5(Cl)Cl'], dtype=object)
test_mols = [Chem.MolFromSmiles(compound)
              for compound in test_dataset.ids]
display_images(mols_to_pngs(test_mols[:5], basename="test"))
1.0s
Python
→

가장 일반적인 scaffolds가 훈련 데이터에 저장되고, 드문 scaffolds는 검증과 테스트 데이터에 저장됩니다.

기계학습 알고리즘의 성능은 데이터의 전처리에 매우 민감할 수 있습니다.

데이터에 적용되는 일반적인 변환 중 하나는 데이터를 평균화하고 단위 표준 편차를 갖도록 정규화하는 것입니다. 이 변환을 로그 용해도에 적용합니다.(-12~2까지의 로그 용해도 분포에 적용)

정규화(normalization)는 다음과 같은 공식을 사용해서 특성 값의 범위를 [0, 1]로 옮깁니다. 이 공식을 사용하면 한 특성 내에 가장 큰 값은 1로, 가장 작은 값은 0으로 변환됩니다. 이 공식을 이용해서 각 특성들의 값을 변환해주면 모두 [0, 1]의 범위를 갖게 됩니다. 이제야 특성들이 평등한 위치에 놓여진 것입니다. 

formula not implementedformula not implemented

반면, 표준화(standardization)는 다음과 같은 공식으로 특성들의 값을 변환해 줍니다.  

여기서 μ는 한 특성의 평균값이고, σ는 표준 편차 입니다. 표준화는 어떤 특성의 값들이 정규분포, 즉 종모양의 분포를 따른다고 가정하고 값들을 0의 평균, 1의 표준편차를 갖도록 변환해주는 것입니다. 표준화를 해주면 정규화처럼 특성값의 범위가 0과 1의 범위로 균일하게 바뀌지는 않습니다.

#데이터 표준화 과정 수행
#출력할 데이터의 범위가 넓기 때문에 데이터의 범위를 조정해줍니다.
transformers = [
  dc.trans.NormalizationTransformer(transform_y=True, dataset=train_dataset)]
for dataset in [train_dataset, valid_dataset, test_dataset]:
  for transformer in transformers:
      dataset = transformer.transform(dataset)
0.9s
Python
→

데이터 처리 후 다음 단계는 간단한 학습 모델을 데이터에 fitting해보는 것입니다. deepchem은 많은 기계 학습 모델 클래스를 제공합니다. 특히, deepChem은 SklearnModel에서 사용할 수 있는 모든 기계 학습 모델을 포함하는 편리한 클래스인 scikit-learn을 제공합니다.

결과적으로 계산된 ECFP4 특징에서 로그 용해도를 예측하는 간단한 랜덤 포레스트 회귀 분석을 시작할 수 있습니다. 모델을 훈련시키기 위해 SklearnModel 객체를 인스턴스화 한 다음 위에서 구축한 train_dataset에서 fit() 메서드를 호출합니다. 그런 다음 모델을 디스크에 저장합니다.

DeepChem 라이브러리를 이용한 머신러닝 모델 학습(RandomForest)

!pip install gast==0.2.2
3.6s
Python
→
from sklearn.ensemble import RandomForestRegressor
sklearn_model = RandomForestRegressor(n_estimators=100)
model = dc.models.SklearnModel(sklearn_model)
model.fit(train_dataset)
4.3s
Python
→

다음으로 검증 세트에서 모델을 평가하여 예측력을 확인합니다. deepchem은 프로세스를 용이하게 하기 위해 평가자 클래스인 Evaluator을 제공합니다. 생성된 모델 객체를 평가하려면, 새로운 Evaluator 인스턴스를 만들고 compute_model_performance()를 호출합니다.

from deepchem.utils.evaluate import Evaluator
metric = dc.metrics.Metric(dc.metrics.r2_score)
evaluator = Evaluator(model, valid_dataset, transformers)
#R2 score는 통계학에서는 결정계수라고 불리는데, 통계학 분야에서 의미하는 바와 같다.
r2score = evaluator.compute_model_performance([metric])
print(r2score)
0.3s
Python
→
from deepchem.utils.evaluate import Evaluator
metric = dc.metrics.Metric(dc.metrics.rms_score)
evaluator = Evaluator(model, valid_dataset, transformers)
rmsScore = evaluator.compute_model_performance([metric])
print(rmsScore)
0.4s
Python
→

기본 랜덤포레스트 회귀 모델은 성능이 매우 좋지 않습니다. 좋은 모델을 구성하기 위해 하이퍼파라미터를 최적화 해봅시다.(모델 사양에서 선택한 항목) 랜덤포레스트에서는 트리 수를 제어하는 n_estimators과 트리 분할을 수행할 때 고려할 특징의 수를 제어하는 max_features를 하이퍼파라미터로 사용합니다. 이제 n_estimators와 max_features를 다르게 선택할 수 있는 리스트를 만들고, 일련의 SklearnModel을 빌드하고 유효성 검사 세트의 성능을 평가해 봅시다.

#dc.models.SklearnModel메서드에 주어진 모델 RandomForestRegressor 및 model_dir를 입력으로 줌
def rf_model_builder(model_params, model_dir):
  sklearn_model = RandomForestRegressor(**model_params)
  return dc.models.SklearnModel(sklearn_model, model_dir)
params_dict = {
    "n_estimators": [10, 100],
    "max_features": ["auto", "sqrt", "log2", None]
}
metric = dc.metrics.Metric(dc.metrics.r2_score)
#dc.hyper.HyperparamOpt를 사용하여 하이퍼 파라미터 최적화
#메서드 안에는 위에서 정의한 rf_model_builder를 입력
optimizer = dc.hyper.HyperparamOpt(rf_model_builder)
#dc.hyper.HyperparamOpt.hyperparam_search메서드에 파라미터 사전, 훈련 데이터세트, 검증 데이터세트, 변환기, 성능측정방법을 입력으로 줌
#최고 성능의 모델, 최고 성능의 하이퍼파라미터, 모든 결과를 변수에 담음
#출력문은 hyperparam_search메서드에서 주어진 log출력물
best_rf, best_rf_hyperparams, all_rf_results = optimizer.hyperparam_search(
    params_dict, train_dataset, valid_dataset, transformers,
    metric=metric)
11.4s
Python
→
#최고 성능을 보이는 모델 하이퍼파라미터 출력
print(best_rf)
#탐색할때 사용한 하이퍼파라미터의 최적 값 출력
print(best_rf_hyperparams)
#모든 탐색 성능 평가 결과 출력
print(all_rf_results)
0.5s
Python
→

최고의 모델은 우리가 구축 한 첫 번째 모델보다 검증 세트에서 높은 inline_formula not implemented를 달성합니다.

하이퍼 파라미터 검색을 수행하지만 간단한 딥 네트워크를 사용합니다.

DeepChem 라이브러리를 이용한 머신러닝 모델 학습(DNN)

!pip3 install 'gast==0.3.2'
3.1s
Python
→
import numpy.random
#random.uniform -자 함수이지만 변동 폭을 조정 low, high, 발생할 갯수
#params_dict호출될때마다 새로운 난수 생성
params_dict = {"learning_rate": [10**-5, 10**-4],
                "decay": [0.01],
               "nb_epoch": [20,50] }
n_features = train_dataset.get_data_shape()[0]
#ResNet Model
# https://arxiv.org/pdf/1603.05027.pdf
def model_builder(model_params, model_dir):
  model = dc.models.MultitaskRegressor(
    1, n_features, layer_sizes=[1000], dropouts=[.25],
    batch_size=50, **model_params)
  return model
optimizer = dc.hyper.HyperparamOpt(model_builder)
best_dnn, best_dnn_hyperparams, all_dnn_results = optimizer.hyperparam_search(
    params_dict, train_dataset, valid_dataset, transformers, 
    metric=metric)
7.7s
Python
→
params_dict
0.0s
Python
→
{'learning_rate': [1e-05, 0.0001], 'decay': [0.01], 'nb_epoch': [20, 50]}
print(best_dnn)
print(best_dnn_hyperparams) 
print(all_dnn_results)
0.3s
Python
→

이제 하이퍼 파라미터를 적절히 선택 했으므로 테스트 세트에서 최상의 모델의 성능을 평가 해 보겠습니다.

rf_test_evaluator = Evaluator(best_rf, test_dataset, transformers)
rf_test_r2score = rf_test_evaluator.compute_model_performance([metric])
print("RF Test set R^2 %f" % (rf_test_r2score["r2_score"]))
0.4s
Python
→
dnn_test_evaluator = Evaluator(best_dnn, test_dataset, transformers)
dnn_test_r2score = dnn_test_evaluator.compute_model_performance([metric])
print("DNN Test set R^2 %f" % (dnn_test_r2score["r2_score"]))
0.3s
Python
→

그래프 컨볼루션 모델 사용해보기

#특성화를 GraphConv모드로 설정하여 데이터 로드
tasks, gcndatasets, gcntransformers =dc.molnet.load_delaney(featurizer='GraphConv')
#데이터 분할
gcn_train_dataset, gcn_valid_dataset, gcn_test_dataset = gcndatasets
#하이퍼파라미터 dict 구성
params_dict = {"n_tasks": [1],
               "mode": ["regression"],
               "dropouts": [0.2],
               "batch_size": [100],
               "nb_epoch": [100,200,300]
               }
def model_builder(model_params, model_dir):
  model = dc.models.GraphConvModel(**model_params)
  
  return model
optimizer = dc.hyper.HyperparamOpt(model_builder)
best_gcn, best_gcn_hyperparams, all_gcn_results = optimizer.hyperparam_search(
    params_dict, gcn_train_dataset, gcn_valid_dataset, gcntransformers, 
    metric=metric)
66.6s
Python
→
print(best_gcn)
print(best_gcn_hyperparams) 
print(all_gcn_results)
0.3s
Python
→
gcn_test_evaluator = Evaluator(best_gcn, gcn_test_dataset, gcntransformers)
gcn_test_r2score = gcn_test_evaluator.compute_model_performance([metric])
gcn_test_r2score
print("GCN Test set R^2 %f" % (gcn_test_r2score["r2_score"]))
0.5s
Python
→
print(gcn_train_dataset.ids[0])
# https://deepchem.readthedocs.io/en/latest/dataclasses.html
print(gcn_train_dataset.X[0])
print(gcn_train_dataset.y[0])
0.7s
Python
→

이제 생성 된 모델에 대해 예측 된 inline_formula not implemented 점수와 실제 inline_formula not implemented 점수를 플로팅하겠습니다.

task = "measured log solubility in mols per litre"
predicted_test = best_rf.predict(test_dataset)
true_test = test_dataset.y
plt.scatter(predicted_test, true_test)
plt.xlabel('Predicted log-solubility in mols/liter')
plt.ylabel('True log-solubility in mols/liter')
plt.title(r'RF- predicted vs. true log-solubilities')
plt.show()
0.5s
Python
→
task = "measured log solubility in mols per litre"
predicted_test = best_dnn.predict(test_dataset)
true_test = test_dataset.y
plt.scatter(predicted_test, true_test)
plt.xlabel('Predicted log-solubility in mols/liter')
plt.ylabel('True log-solubility in mols/liter')
plt.title(r'DNN predicted vs. true log-solubilities')
plt.show()
1.1s
Python
→
predicted_test = best_gcn.predict(gcn_test_dataset)
true_test = gcn_test_dataset.y
plt.scatter(predicted_test, true_test)
plt.xlabel('Predicted log-solubility in mols/liter')
plt.ylabel('True log-solubility in mols/liter')
plt.title(r'GCN predicted vs. true log-solubilities')
plt.show()
0.6s
Python
→
Runtimes (1)