Attention-based Multi-Input Deep Learning Architecture for Biological Activity Prediction: An Application in EGFR Inhibitors

TL;DR

This paper introduces an attention-based deep learning model that simultaneously uses structural and chemical data to predict biological activity of EGFR inhibitors, improving accuracy and interpretability.

Contribution

The study presents a novel multi-input architecture with attention mechanism that enhances prediction performance and interpretability in drug activity modeling.

Findings

Achieved MCC of 0.58 and AUC of 90% on EGFR dataset.

Outperformed existing models in bioactivity prediction.

Enabled interpretation of chemical structure contributions.

Abstract

Machine learning and deep learning have gained popularity and achieved immense success in Drug discovery in recent decades. Historically, machine learning and deep learning models were trained on either structural data or chemical properties by separated model. In this study, we proposed an architecture training simultaneously both type of data in order to improve the overall performance. Given the molecular structure in the form of SMILES notation and their label, we generated the SMILES-based feature matrix and molecular descriptors. These data were trained on a deep learning model which was also integrated with the Attention mechanism to facilitate training and interpreting. Experiments showed that our model could raise the performance of prediction comparing to the reference. With the maximum MCC 0.58 and AUC 90% by cross-validation on EGFR inhibitors dataset, our architecture was outperforming the referring model. We also successfully integrated Attention mechanism into our model, which helped to interpret the contribution of chemical structures on bioactivity.

Figures

Tables

Table 1. TABLE I: Features table

Features	Description	Size
Type of atom	H, C, O, N, or others	5
No. of Hs	Total number of attached Hydrogen atoms	1
Degree	Degree of unsaturation	1
Charge	Formal charge	1
Valence	Total valence	1
Ring	Included in a ring or not?	1
Aromaticity	Included in a aromatic ring or not?	1
Chirality	R, S, or others	3
Hybridization	$s,sp,s{p}^2,s{p}^3,s{p}^{3}d,s{p}^3{d}^2,$ or others	7
SMILES symbol	( ) [ ] . : = # \ / @ + -	21
	Ion_charge Start End
Total		42

Table 2. TABLE II: Hyper-parameter in CNN branch

Layer	Hyper-parameter	Value
$1^{st}$ conv2d	No. of input channels	1
	No. of output channels	6
	Kernel size	(3,3)
$2^{nd}$ conv2d	No. of input channels	6
	No. of output channels	16
	Kernel size	(3,3)
Dropout	Dropout rate	Be tuned

Table 3. TABLE III: Hyper-parameter in MD branch

Layer	Hyper-parameter	Value
$1^{st}$ Linear layer	No. of neurons	512
$2^{nd}$ Linear layer	No. of neurons	128
$3^{rd}$ Linear layer	No. of neurons	64

Table 4. TABLE IV: Hyper-parameter tuning

Hyper-parameter	Value	Optimal value
Batch size	32; 64; 128; 512	128
Dropout rate	0; 0.2; 0.5	0.5
Optimizer	SGD; ADAM	ADAM
Learning rate	1e-4; 1e-5; 1e-6	1e-5
Threshold	0.2; 0.5; 0.8	0.2

Table 5. TABLE V: Performance Metrics

Metrics	Formulas†
Sensitivity	$\frac{\text{TP}}{\text{TP}+\text{FN}}$
Specificity	$\frac{\text{TN}}{\text{TN}+\text{FP}}$
Accuracy	$\frac{\text{TP}+\text{TN}}{\text{TP}+\text{TN}+\text{FP}+\text{FN}}$
MCC	$\frac{\text{TP}\times \text{TN}-\text{FP}\times \text{FN}}{\sqrt{\left(\text{TP}+\text{FP}\right)\left(\text{TP}+\text{FN}\right)\left(\text{TN}+\text{FP}\right)\left(\text{TN}+\text{FN}\right)}}$
AUC	the area under the ROC curve

Table 6. TABLE VI: Performance Comparison

Metrics †	H. Singh	CNN	CNN + MD	CNN + MD
	et al.			+ ATT ‡
LOSS	N.A	0.2727	0.2434	0.2399
SENS	69.89%	74.31%	75.29%	74.11%
SPEC	86.03%	85.77%	89.75%	90.15%
ACC	83.66%	84.10%	87.66%	87.83%
MCC	0.49	0.50	0.58	0.57
AUC	89.00%	87.52%	90.32%	90.84%

Equations

$$y=f\left(x_1,x_2,...,x_n;\theta \right)+ϵ$$

$$\hat{y}=\sigma \left(w_0+x_1{w}_1+x_2{w}_2+...+x_n{w}_n\right)$$

$$l_i^{(1)}=\sigma \left(w_{10}^{(1)}+\sum w_{1i}^{(1)}{x}_i\right)$$

$$l_i^{(2)}=\sigma \left(w_{20}^{(2)}+\sum w_{2i}^{(2)}{l}_i^{(1)}\right)$$

$$\hat{y}=\sigma \left(w_{30}^{(3)}+\sum w_{3i}^{(3)}{l}_i^{(2)}\right)$$

$$O_{i,j}^k=\sum _{i^{\prime }=1}^w\sum _{j^{\prime }=1}^w\sum _{l=1}^c[X{]}_{i+i^{\prime }-1,j+j^{\prime }-1,l}[F_k{]}_{i^{\prime },j^{\prime },l}$$

$$\hat{\theta }=\arg \underset{\theta }{\min }J(\theta )=\frac{1}{n}\sum _{i=1}^n\mathcal{L}(x_i,y_i;\theta )$$

$${\theta }_{t+1}={\theta }_t-\gamma {\nabla }_{\theta }J(\theta )$$

$${\hat{x}}_i=\frac{x_i-{\mu }_{\mathcal{B}}}{\sqrt{{\sigma }_{\mathcal{B}}^2+ϵ}}$$

$${\hat{x}}_i=\gamma {\hat{x}}_i+\beta$$

$${\vec{a}}_i=\frac{1}{1+\exp \left({\vec{R}}_i\cdot \vec{m}\right)}$$

$$\vec{f}=\sum {\vec{a}}_i\times {\vec{R}}_i$$

Full text

††publicationid: pubid: 978-1-7281-3003-3/19/$31.00 ©2019 IEEE

I Introduction

Machine learning was applied widely in drug discovery, especially in virtual screening for hit identification. The most popular techniques are Support Vector Machine (SVM), Decision Tree (DT), k-Nearest Neighbor (k-NN), Naive Bayesian method (NB), and Artificial Neural network (ANN). [1]. In these methods, ANNs need not assume that there was any type of relationship between activity and molecular descriptors, and ANNs are usually outperforming in traditional Quantitative structure – activity relationship problem because they can deal with both nonlinear and linear relationship. As a result, ANNs rose to become a robust tool for Drug Discovery and Development [2]. However, ANNs are usually sensitive to overfitting and difficult to design an optimal model. Additionally, ANNs also require huge computation resources and their results usually are unable to be interpreted. Those weaknesses could be a reason for limited use of neural network comparing to Decision Tree or Naive Bayesian algorithms [2, 1].

Above algorithms can be applied to various types of chemical features. These features are either structural information or chemical properties. The structural information could be represented as fingerprint vector by using specific algorithms (e.g, Extended-Connectivity Fingerprints [3], Chemical Hashed Fingerprint [4]) while chemistry information could be described by various molecular descriptors (e.g, logP, dipole moment). The ideas that combine some types of features to improve the overall performance was also mentioned in a number of studies [5, 6]. In these models, each set of chemical features was trained by specific algorithms (SVM, DT, k-NN, NB, ANN or any other algorithms) to generate one particular output. After that, these outputs were pushed to the second model which was usually another multi-layers perceptron model before giving the final result. The disadvantage of this approach is that we need to train each feature set separately because of completely different algorithms. As a result, it is difficult to build a pipeline for all training algorithms. In other words, the automation of the training procedure was reduced.

Regarding the interpretation of neural network model, there are some interests in making neural network models more explainable and interpretable. A worthy approach needs to be mentioned is the use of Attention mechanisms in sequence-to-sequence model. With the encoder-decoder architecture, the attention approach not only improves the performance but reveals the alignment between input and output [7].

To deal with the problem of both automation and interpretation in predicting biological activity, we made an effort to combine different types of chemical features in one deep learning architecture and also integrate Attention mechanism. As a result, our model could train concurrently several feature sets and explain the interaction between the features and outcomes.111Our models are available at https://github.com/lehgtrung/egfr-att

II Background

II-A Overview of Neural network

II-A1 Artificial neural network

The Artificial neural network is computing architecture which enables a computer to learn from historical data. Nowadays, it is one of the main tools used in machine learning. As the name suggests, artificial neural networks are inspired by how biological neurons work, however, an artificial neural network is a composition of many differentiable functions chained together. Mathematically, a neural network is a non-linear mapping which assumes the output variable $y$ as a non-linear function of its input variables $x_{1},x_{2},...x_{n}$

[TABLE]

where $\theta$ is the parameters of the neural network and $\epsilon$ is model’s inreducible error.

A very simple neural network which contains only input and output is described as follows:

[TABLE]

where $\hat{y}$ is an approximation of $y$.

As shown in Fig. 1, each input variable $x_{i}$ is represented as a node in the input layer and connects to the output node through a connection with weight $w_{i}$. Note that a node with value $1$ is attached to the input layer with the corresponding connection weight $w_{0}$ to represent the offset term. The $\sigma$ function is called activation function or squash function which introduce non-linear relationship between input $\mathbf{x}$ and output $y$. The Sigmoid $\left(f(x)=\frac{1}{1+e^{-x}}\right)$ and ReLU $\left(f(x)=max(0,x)\right)$ function are the most widely used activation functions. The model in Fig. 1 is referred as generalized linear model.

The generalized linear model is simple, thus may not be able to describe complex relationship between inputs and outputs. Therefore, this architecture can be extended by building multiple generalized linear model as the form of layers (or fully connected layers or hidden layers) and stacking those layers together to build a neural network. Fig. 2 illustrates a two layered neural network, we also add neuron $1$ to the second layer as we do for the generalized linear model below. Let $l^{(k)}$ and $\hat{y}$ represent the $k$-th hidden layer and output layer respectively, neural network in Fig. 2 can be described mathematically as follows:

[TABLE]

[TABLE]

[TABLE]

Note that a neural network may consist of an arbitrary number of many layers by simply stacking more layers. A network containing more than one layers is usually called a deep neural network.

II-A2 Convolutional neural network

Convolutional neural network (CNN) [8] is a class of neural network. The models using CNN are usually designed to operate on data with grid-like topology. CNN models are usually considered as the state-of-the-art architectures in the computer vision related tasks. CNNs are also applied in biological tasks and achieved remarkable results [9, 10].

Basically, a CNN block is a combination of convolution layers followed by non-linear activation and pooling layers.

• •

Convolutional layer (CONV): A convolutional layer is composed of a set of kernels to extract local features from the previous input. Each kernel is represented as a 3D tensor $F_{k}\in\mathbb{R}^{w\times w\times c}$, where $w$ is the size of the kernel (typically $3$ or $5$) and $c$ denotes the total number of kernels. Since $c$ is equal to the input’s third dimension, it is frequently omitted when referring to the kernel shape. For an input $X\in R^{d_{1}\times d_{2}\times d_{3}}$, each kernel $F_{k}$ convolves with $X$ to attain a single feature map $O^{k}\in\mathbb{R}^{((d1-w+)/s+1)\times(d2-w)/s+1)}$ where

[TABLE]

where $[X]_{i,j}\in\mathbb{R}^{w\times w\times d_{3}}$ is a small block (known as receptive field) of $X$ around location $(i,j)$; $s$ is the stride which is the interval of the receptive fields of neighboring units.

• •

Pooling layer (POOL): Pooling layer creates a summary of learned features from CONV layers by aggregating the information of nearby features into a single one. The most common design of pooling operation is max-pooling. For example, a $2\times 2$ max-pooling filter operating on a particular feature map $F$ with size $(n,m)$ will compute $\max\{F_{i,j},F_{i+1,j},F_{i,j+1},F_{i+1,j+1}\}$ for each coordinate $i,j$ in $F$. This will result in a new features map with size $(m/2,n/2)$. Since a CNN typically contains multiple stacking of CONV layers, pooling is used to reduce data dimension which causes the model less computationally expensive. Pooling can also make the model invariant to small positional and translational changes.

A typical CNN architecture is generally made up of series of CNN blocks followed by one or more fully connected layers at the end. Fig. 3 illustrates a simple CNN architecture for image classification problem.

II-A3 Training Neural Networks

The goal of learning is to minimize the loss function with respect to the network parameters $\theta$. To do that, we need to find an estimate for the parameters $\hat{\theta}$ by solving an optimization problem of the form

[TABLE]

where $n$ is the number of instances in the training dataset; $\mathcal{L}$ is the loss function which measures the discrepancy between model output and the ground truth. Because the optimization problem does not have closed form solution, the method of gradient descent is used. Firstly, the parameters $\theta$ are randomly constructed, for every iteration, the parameters are updated as follow

[TABLE]

This process continues until some criterion is satisfied. Here, $\gamma$ is a constant called the learning rate which is the amount that the weights are updated during training. As presented in the equation 4, the loss function is computed over all the examples, which is computationally extensive. In practice, we use a modified version of gradient descent called stochastic gradient descent , that means, we do not use the whole dataset for gradient computation but a subset of data called a mini-batch. Typically, a mini-batch contains from dozens to hundreds of samples depending on system memory. Since the neural network is a composition of multiple layers, the gradient with respect to all the parameters can be methodically computed using the chain rule of differentiation also known as back-propagation algorithm.

II-A4 Regularization

One of the major issues in training neural networks is overfitting. Overfitting happens when a network performs too well on the data it has been trained on but poorly on the test set which it has never seen before. This phenomenon is due to the large number of parameters in the network. Regularization is able to regulate a network activity to ensure the model actually learns the underlying mapping function not just memorizing the input and output. Recently, there are two advanced regularizers which are widely used in the deep neural network.

• •

Dropout: During training, some weights in the network at a particular layer could be co-adapted together which may lead to overfitting. Dropout tackles this issue by randomly skipped some weights (explicitly set them zero) with a probability $p$ (usually $p=0.5$ or $0.8$). During inference, dropout is disabled and the weights are scaled with a factor of $p$ [11].

• •

Batch normalization: Recall in regression analysis, one often standardizes the designed matrix so that the features have zero mean and unit variance. This action called normalization speeds up the convergence and make initialization easier. Batch normalization spread this procedure to not only input layer but all of the hidden layers. During training, let $x_{i}$ is values across a mini-batch $\mathcal{B}=\{x_{1},x_{2},...,x_{k}\}$, the batch norm layer calculate normalized version $\hat{x}_{i}$ of $x_{i}$ via:

[TABLE]

where $\mu_{\mathcal{B}}=\frac{1}{k}\sum_{i=1}^{k}x_{i}$; $\sigma^{2}_{\mathcal{B}}=\frac{1}{k}\sum_{i=1}^{k}(x_{i}-\mu_{\mathcal{B}})^{2}$ are mini-batch mean and variance respectively, $\epsilon$ is a constant to help computational efficiency. To make it more versatile, a batch norm layer usually has two additional learnable parameters $\gamma$ and $\beta$ which stand for scale and shift factor such that:

[TABLE]

During inference, mini-batch mean and variance are replaced by population mean and variance which are estimated during training [12].

II-A5 Attention mechanism

Neural networks could be considered as a ”black box” optimization algorithm since we do not know what happens inside them. Attention mechanism enables us to visualize and interpret the activity of neural networks by allowing the network to look back to what it has passed through. This mechanism is motivated by how we, human, pay visual attention to certain regions of images or important words in a sentence. In the neural network, we can simulate this behavior by putting attention weights to express the importance of an element such as pixel in an image or a word in a sentence.

Attention mechanism was applied widely and now becomes a standard in many tasks such as Neural machine translation [13], Image captioning [14].

II-B Overview of EGFR

Epidermal Growth Factor Receptor (EGFR) is a member of ErbB receptor family that consists of 4 types: EGFR, HER1, HER2/new, HER3, and HER4. They are located in the cell membrane with the intrinsic tyrosine kinase. The binding of ligands (TGF-$\alpha$, amphiregulin, and other ligands) and EGFR triggers the signal amplification and diversification which lead to cell proliferation, apoptosis, tumor cell mobility, and angiogenesis. In some type of cancer (such as lung cancer), the overexpression and constitute activation cause the dysregulation of EGFR pathway that activates the tumor process [15, 16].

For two decades, there was a great deal of effort in studying this target to discover novel medicine. 3D-QSAR was studied widely to analysis the molecular filed of ligands, which reveals the relationship between various substituents on molecules and biological activity [17, 18, 19, 20, 21, 22]. Other methods were also useful. R. Bathini et al. employed the molecular docking and molecular mechanics with generalized born surface area (MM/GBSA) to calculate the binding affinities of protein-ligand complexes [19]. G. Verma et al. conducted pharmacophore modeling in addition to 3D-QSAR to generate a new model which was used for screening novel inhibitors [22].

Regarding the application of machine learning techniques in EGFR inhibitors discovery, H. Singh et al. [23] used Random Forest algorithms to classify EGFR inhibitors and non-inhibitors. In their study, the authors collected a set of diverse chemical and their activity on EGFR. A model with high accuracy was trained and validated by 5-fold cross-validation (0.49 in MCC and 0.89 in AUC).

II-C Overview of Features set

II-C1 SMILES Feature matrix

SMILES (Simplified Molecular Input Line Entry System) is a way to represent the molecules in in silico study. This method uses a strict and detailed set of rules to interpret the molecular structure into the chemical notation which is user-friendly but also machine-friendly. In particularly, SMILES notation of a molecule is a chain of character which is specified for atoms, bonds, branches, cyclic structures, disconnected structures, and aromaticity [24].

Based on this representation, M. Hirohara et al. developed a SMILES-based feature matrix to train a convolutional neural network model for predicting toxicity. His model outperformed the conventional approach and performed comparably against the winner of Tox21 data challenge [25]. In this dataset, each molecule was represented in the form of SMILES notation and the output consisted of 12 tasks to predict [26].

II-C2 Molecular descriptors

Molecular descriptors are terms that characterize a specific aspect of a molecular, including substituent constants and whole molecular descriptors [27]. The calculation of former type derived from the difference in functional group substitution into the main core of the compound. Based on this approach, the latter is the expansion of the substituent constant. However, some whole molecular descriptors are developed from totally new methods or based on physical experiments [28].

III Method

III-A Dataset

The dataset used in this study was collected by H. Singh et al. [23]. This dataset contains 3492 compounds which is classified as inhibitor or non-inhibitor of EGFR. The inhibition activity of a particular substance was assigned if its $IC50$ is less than 10 $nM$. The ratio of inhibitors over non-inhibitor is $506:2986\approx 1:6$. The information of chemical includes ID, SMILES representation, and class (1 for inhibitor and 0 for non-inhibitor).

III-A1 SMILES Feature matrix generation

Based on the collected dataset, the chemical structure data in the form of SMILES notation was preprocessed and converted to the canonical form which is unique for each molecule by the package rdkit [29].

In this study, the SMILES Feature matrix generation method developed by M. Hirohasa et al. was used to encode the chemical notation. The maximum length of each input was 150 and thus the input strings with length below 150 were padded with zeros at the tail. In their method, for each character in the SMILES string, a 42-dimensional vector was computed (Table I). The first 21 features represent the data about the atom and the last 21 features contain SMILES syntactic information [25].

III-A2 Descriptor calculation

We used the package mordred built by H. Moriwaki, Y. Tian, N. Kawashita et al. [30] to generate molecular descriptor data. Because the SMILES notation do not provide exact 3D conformation, the 2D descriptors were only calculated with total of 1613 features . The generated data was preprocessed by imputing the meaningless features or the variables which are same for whole dataset. A standard scaler was also used to normalize the molecular descriptors dataset. We used package numpy [31], pandas [32] and scikit-learn [33] for this process.

III-B Model architecture

III-B1 Convolutional neural network (CNN) branch

The SMILES Feature matrix was flown through 2 CNN blocks each consisted of a 2D convolution layer, one Batch normalization layer, a Dropout layer, and a Max pooling layer before being flattened and fully connected via a linear layer. The detail of the hyper-parameters of each layer is represented in Table II.

III-B2 Attention mechanism

The idea of using the attention mechanism came from the fact that each chemical’s atoms contribute differently to the drug’s effect. In other words, we put attention or weight to the atoms in the chemical which are represented by rows in the SMILES feature matrix. The larger the weight of an atom is, the more contribution of which atom contribute to the drug. By doing this, we can extract each components weights for interpreting the results and analysis.

Let $\vec{m}$ and $R$ denote the vector obtained by feeding SMILE vector through a linear layer and the SMILE feature matrix, respectively. The model uses the similarity between $\vec{m}$ and the $i$-th row $\vec{R}_{i}$ as a measure of the importance of $\vec{R}_{i}$. Particularly, let $\vec{a}_{i}$ denote the attention weight vector, it is formulated as follows:

[TABLE]

After that, the vector $\vec{a}$ is then used as the coefficient of a linear combination of rows in the SMILES feature matrix. Therefore, the output of the attention layer $\vec{f}$ is expressed as:

[TABLE]

III-B3 Molecular Descriptors (MD) branch

The MD branch is used to train the molecular descriptors data. Let $\vec{t}$ is the vector obtained by feeding molecular descriptors vector through 3 blocks of fully connected layers each consisted of a fully connected layer, a batch normalization layer, and a dropout layer. $\vec{t}$ is considered as a high-level representation of molecular descriptors data.

III-B4 Concatenation

The three vectors: $\vec{f}$, $\vec{m}$, $\vec{t}$ are then concatenated as in Fig. 5 to form the final representation vector. This vector combines information from both molecule structural information and descriptors information. It is then fed through a linear layer and squash by sigmoid function to make the final prediction as the probability of the molecule as an inhibitor.

III-C Hyper-parameter tuning

In this study, PyTorch platform [34] was used in order to implement our model and the 5-fold cross-validation was conducted to evaluate the performance. In this method, the dataset was split into 5 parts. For each fold, the model was trained on the set of 4 parts and tested on the remaining part. The choice of the training set was permuted through all divided parts of the dataset, thus the model was trained 5 times and the average performance metrics of each time was used to evaluate. The ending point of Training step was determined by Early-stopping technique [35, 36]. Thus, for each fold, the model would stop training if the loss value increases continuously 30 epochs.

The second column of Table IV presents the hyper-parameters and their considered values in the tuning step. Grid Search technique was conducted to determine the best combination of hyperparameters which had the best performance. However, in the case of discovering the suitable batch size for training, several suggested values were tested and the chosen was the value which utilized the system efficiently. Additionally, the threshold of the classifier was determined by analyzing the ROC plot and the Precision-Recall curve. The most optimal threshold was the point nearest to top-left of the ROC plot and gave the balance between Precision and Recall in the latter plot.

III-D Performance Evaluation

In order to assess the performance of each model, several metrics were calculated during training and validation steps (Table V).

The ROC analysis and AUC are usually considered as the most popular metrics for imbalanced dataset because they are not biased against the minor label [37, 38]. However, these metrics show the overall performance in the whole domain of threshold. In other words, ROC or AUC are not represented for a particular classifier. In our study, MCC was the most preferred criteria to evaluate model performance. This is because MCC considers all classes in the confusion matrix whereas other metrics (eg, accuracy or F1-score) do not fully use four classes in the confusion matrix [39]. The remaining metrics were still useful for the benchmark.

IV Result

IV-A Hyper-parameters optimization

Despite of imbalanced classes, the loss function of model (binary cross-entropy) still converged with the loss value of $0.10-0.16$ for Training set and $0.22-0.28$ for Validation set at the end of each fold when cross-validating the model CNN + MD + ATT

The training batch size is 128 which gave the best utility on the GPU Tesla K80. When comparing the effect of two types of the optimizer, we observed that Adaptive Moment Estimation (ADAM) showed a better result than Stochastic gradient descent (SGD) in both running time and model performance.

The optimal collection of hyper-parameters was listed in the third column of Table IV. These values were used for evaluating the performance of three considered models.

IV-B Performance

Our architecture was trained on the EGFR dataset and evaluated by mentioned cross-validation method. The final result represents in Table VI. When using structure information in the CNN branch only, the performance was slightly better than H. Singh et al. model in all metrics excepting AUC. However, the CNN + MD model and CNN + MD + ATT model was outperforming to the reference model. By comparing two important indicators, it can be seen that there was an improvement in MCC and AUC. Additionally, when training with more branch (CNN + MD and CNN + MD + ATT), the Running time was also reduced significantly to around a half.

IV-C Attention mechanism

The Attention mechanism was successfully implemented in our architecture. From the attention weights vector, the weights representing to each atom in the molecules were extracted and used to indicate their contribution to the compound’s activity. By using visualization on the package rdkit, these attention weights could be used to visualize the distribution of contribution over the molecular structure. Figure 7 illustrates some example from the model.

Additionally, the integration of attention mechanism to the model did not give much more improvement in performance. We believed that this problem could be improved by optimizing the hyperparameters. However, we did not focus on hyperparameters optimization for this model because the aim of Attention vector was to make the model more interpretable.

The visualization of Attention vectors revealed some key findings in chemical structure, such as the Nitrogens in the hetero-cyclic structure are usually highlighted by the model and the halobenzyl substitution on heterocyclic structure contributed positively to the bioactivity.

V Discussion

There are two major advantages to our architecture. The first strong point is the combination of both structure information and chemical attributes in a single learning model. As a result, this advancement made a significant improvement in both performance and automation. Another worthy innovation was the integration of attention mechanism which facilitated the interpretation of the model. In fact, attention weight generated by the model would help explain the contribution of each atom on the overall biological activity.

Comparing to another effort to make deep learning model more interpretable, our architecture has an advantage in computation because it is easier to generate the SMILES Feature matrix than other algorithms. For example, Sanjoy Dey et al. [40] used ECFP fingerprint algorithms to transform the chemical structure into the matrix feature. This method does not treat the molecule as a whole structure but calculate on each fragment of a chemical with a particular radius. Additionally, there are required calculation to generate the features including tuning the hyper-parameters of the algorithms (e.g the radius of calculation).

Regarding to our implementation of Attention mechanism, we observed that each atom in a substance was treated separately; as a result, the connection between atom was not highlighted in our model, as well as the contribution of some functional groups which contain many atoms (e.g, carbonyl, carboxylic, etc) was not clearly illustrated. We proposed a solution for this limitation that is to add another branch to the architecture which embeds the substructure patterns (e.g, Extended-Connectivity Fingerprints [3], Chemical Hashed Fingerprint [4]).

Additionally, the lower performance of the CNN model comparing to the baseline model was another interesting finding. This could be due to the sparsity of SMILES feature matrix. In fact, the CNN as well as other deep learning algorithms require much data to accumulate the information in the training step. In case of SMILES feature data, because of zero paddings to justify the length of encoding vectors, the feature matrix became sparse. This led to the fact that the model required more data for training but the dataset was quite small for deep learning. However, in the case of CNN + MD or CNN + MD + ATT model, because of the addition of another input data, the models acquired information more easily. As a result, the performance was improved in terms of all metrics.

When considering the running time between different models, it is clear that the longest running time was that of model with only CNN branch (59 min) while the more complicated model with more data like CNN + MD and CNN + MD + ATT took just a half of running time with 31 min and 37 min, respectively. This could be because the SMILES Feature matrix in the CNN model was sparse so the model should train longer to achieve the convergence of loss function. In the CNN + MD and CNN + MD + ATT model, there could be a complement between different input branches and we supposed that there was an information flow transferring between two branches, which facilitated the training stage and performance improvement. In other studies which also used several types of data [5, 6], the model trained model separately and did not use this information connection. This phenomenon might represent an advantage of our architecture.

In conclusion, the combination of different source of features is definitely useful for bioactivity prediction, especially when using deep learning model. The attention-based multi-input architecture we proposed achieved a superior score comparing to referring model. Additionally, the attention mechanism would help to interpret the interaction between each element of chemical structures and their activity.