The ATP dataset utilized in this study was constructed by our group (Hu et al., 2020) through the following steps: initially, 1728 ATP protein chains were sourced from the semi-manual BioLip database (Yang et al., 2013); subsequently, these chains were filtered with sequence length (>50 residues), the resolution (<3 Å), and the sequence identity (<30%); resulting in the ATP-289 dataset with 289 protein chains including 3901 ATP binding residues and 104153 ATP non-binding residues; the dataset was then partitioned randomly into training and testing sets, with the former containing 260 protein chains encompassing 3526 ATP binding residues and 92804 non-binding residues, and the latter comprising 29 protein chains with 375 ATP binding residues and 11349 ATP non-binding residues. We performed five-fold cross-validation using the training set to obtain the trained model and validated the model’s effectiveness using an independent test set. The source codes and datasets in this study are available at https://github.com/tlhsx/S-DCNN.

Previous researches indicated that residues neighboring binding sites can influence ligand interactions with these sites (Hu et al., 2016; Liu et al., 2019). To address this, protein sequences were segmented into fragments using the sliding window method, ensuring each amino acid resided at the fragment center by adding (L-1)/2 pseudo-amino acids at both sequence ends. Here, L denotes the fragment length. If a residue of (L + 1)/2 was the binding residue, it was defined as a positive sample, otherwise, it was a negative sample. Based on the previous references (Yu et al., 2013; Ding et al., 2017; Zhao et al., 2019; Hu et al., 2020; Hu et al., 2021), the intercepted fragment L was 17.

The secondary structure of proteins can reflect the trend of the backbone chain, and the dihedral angle is the main descriptor of the secondary structure, which can reflect the local structural information of proteins and is a very effective feature for predicting protein-ligand binding residues (Chen et al., 2011; Cui et al., 2019; Liu et al., 2020). Here, we applied firstly the reclassified dihedral angles to the prediction of ATP binding residues. First, the values of phi (φ) and psi (ψ) angles were obtained from the primary sequence using ANGLOR software (Wu and Zhang, 2008), and the value range of φ and ψ angles both were [−180°, 180°]; then every 15° was divided into an interval, the φ and ψ angles both were divided into 24 intervals; the difference value of probability of the φ and ψ angles between the positive and negative samples was obtained, as shown in Figures 1A, B. The formula for calculating the probability difference is expressed in Equation 1 as follows: Δ P i = P i + − P i − (1) where, P i j = n i j ∑ i = 1 24 n i j , n _ij represents the number of i ^th interval in positive or negative samples; i (i = 1, 2, … , 24) represents the divided interval; j (+ or -) represents a positive or negative sample.

In Figure 1, it was found that there exist significant differences in the probabilities of the φ and ψ angles between the positive and negative samples. Using 0 as the threshold, we divided the φ and ψ angles into three intervals, which were represented by functions g(x) and h(x) (i.e., Equations 2, 3). Through the above analysis, we selected the reclassification information of the φ and ψ angles as features. g x = I , x ∈ − 180 ∘ , − 90 ∘ II , x ∈ − 90 ∘ , − 60 ∘ III , x ∈ − 60 ∘ , 180 ∘ (2) h x = I , x ∈ − 180 ∘ , − 60 ∘ II , x ∈ − 60 ∘ , 0 ∘ III , x ∈ 0 ∘ , 180 ∘ (3)

In accordance with the principles of physics, the stability of molecular structures increases as energy decreases (Wang et al., 2021). In consideration of the specificity of ATP binding to proteins, we analyzed the Laplace energy values of the 20 amino acids between positive and negative samples in Figure 2. The analysis revealed varying energy probabilities among the 20 amino acids between the positive and negative samples. Consequently, the amino acids were regrouped into four categories: the first group comprised G, I, S, T and V, with markedly higher values in the positive set than in the negative set; the second group included C, H and M, where the positive set’s values slightly surpassed those of the negative set; the third group encompassed F, K, N, R, W and Y, in which the values of negative set were slightly higher than that of positive set; the fourth group consisted of A, D, E, L, P and Q, with notably higher values in the negative set than in the positive set. Subsequently, the energy reclassification details were utilized as feature parameters for ATP binding residue identification.

Researchers have analyzed the influence of binding residues and their neighboring residues on the protein-ATP binding process at the sequence fragment level. In protein-ligand interactions, the amino acids’ specific preferences in crucial binding residues that directly engage with the ligand play a vital role in the binding process. Hence, we introduced a novel parameter extraction method termed propensity factors. Originally suggested by Chou and Fasman (1974), propensity factors have found utility in predicting protein secondary structures and ion ligand-binding sites (Chou and Fasman, 1979; Xu et al., 2022). The formula was expressed in Equation 4 as follows: F i j = p i j p j (4) where, p i j = n i j N i , p j = N j N t , N i = ∑ i = 1 20 n i j , N t = ∑ j = 1 2 N j , n _ij represents the number of amino acid i in binding residues or non-binding residues; N _j represents the number of binding residues or non-binding residues; i (i = 1, 2, … , 20) represents 20 amino acids; j (j = 1, 2) represents binding residues and non-binding residues. The propensity factors of 20 amino acids were statistically analyzed, as shown in Figure 3.

The propensity factor values of amino acids D, G, H, K, R, S, and T within the binding residues exhibited notably higher values compared to those in the non-binding residues. Hence, the propensity factors, serving as a novel extraction method, effectively capture the preferences of the binding residues.

Utilizing sequence information, we extracted amino acids and derived predicted secondary structure information, relative solvent accessibility, and the hydrophilic-hydrophobic profile as fundamental features. These parameters, extensively employed in prior research, have demonstrated exceptional predictive capabilities (Chen et al., 2011; Zhang et al., 2012; Yu et al., 2013; Hu et al., 2016; Ding et al., 2017; Zhao et al., 2019; Hu et al., 2020; Song et al., 2020a; Hu et al., 2021). The hydrophilic-hydrophobic properties were used to classify the 20 amino acids into six distinct categories (Pánek et al., 2005). Secondary structure and solvent accessibility predictions were generated through ANGLOR software, categorizing secondary structure into α-helix, β-sheet and coil. Following guidelines from a source (Hu et al., 2020), relative solvent accessibility predictions were partitioned into four intervals: (0, 0.2], (0.2, 0.45], (0.45, 0.6], (0.6, 0.85].

The researchers observed significant disparities in amino acid frequencies between positive and negative samples, prompting the utilization of amino acid composition as a parameter (Hu et al., 2020; Wang et al., 2021; Sun et al., 2022). Here, from the amino acids composition, secondary structure, relative solvent accessibility, φ angle, ψ angle and energy, we extracted 21, 4, 5, 4, 4 and 5-dimensional composition information, respectively.

Previous studies have shown that the position weight matrix can well reflect the site conservation of amino acids in protein sequences (Hu et al., 2020; Xu et al., 2022). Here, the matrix elements were expressed in Equation 5 as follows: m i , j = ln p i , j p 0 , j (5) where, p i , j = n i , j + N i q N i + N i , N i = ∑ j = 1 21 n i , j , P _0,j represents the background probability, and n _i,j represents the frequency of the j ^th amino acid at the i ^th site, j represents 20 kinds of amino acids and vacancies, q represents the number of classifications, here it is 21. Two standard scoring matrices were obtained from the positive and negative training sets, and 2L-dimensional feature vector were obtained for each segment. Similarly, the predicted secondary structure (q = 4), relative solvent accessibility (q = 5), energy (q = 5), φ angle (q = 4) and ψ angle (q = 4) were also extracted by the same methods, and a total of 6 × 2L-dimensional site conservative information was obtained.

The intermolecular hydrophobic effect was a complex process which was mainly determined by the entropy effect (Wang et al., 2021). The information entropy was an effective method to extract information of hydrophilic-hydrophobic (Wang et al., 2021; Sun et al., 2022; Xu et al., 2022). Here, the information entropy formula was expressed in Equation 6 as follows: H x = − ∑ j = 1 q p j ⁡ log 2 ⁡ p j (6) where, p j = n j N , n _j represents the frequency of occurrence of the j ^th classification in a segment, N is the segment length, and q represents the hydrophilic-hydrophobic classifications and vacancies, here it is q = 7.

As one of the most important branches of deep learning framework, DCNN usually consisted of the input layer, convolutional layer, pooling layer, fully connected layer and output layer. The potential complex information was detected for the input raw data, and then through a series of high-dimensional and high-level projection mapping, the deeper representation information of the classified objects was obtained. The alternating distribution of the convolutional layer and the pooling layer made the convolutional neural network had better fault tolerance and parallel processing ability, and the generalization ability and adaptability of the model were greatly enhanced. It has been widely used in various fields.

The DCNN model framework in this paper was implemented by Keras, and the bottom layer was based on the TensorFlow framework. Here, the batch normalization was employed to avoid vanishing gradients and speed up the convergence of the network. In order to prevent over-fitting of the model, the layer of dropout was employed. The relu nonlinear activation function was used to improve the expressive ability of the model and greatly shortened the learning cycle. To effectively avoid over-fitting problems caused by continued training, the early stopping module was used. Adam and cross-entropy were used as optimizer and loss function, respectively. The output layer applied the sigmoid function to make the classification objects output probability values between 0 and 1. The hyperparameters in DCNN algorithm had an influence on the training speed and performance of the predictor. Based on the previous research, we mainly optimized the following three hyperparameters. Here, the dropout was set as 0.2; the range of the number of convolutional layers was from 1 to 6; the range of filters and batch size was both from 2 to 128. Detailed description of DCNN architecture can be viewed at https://github.com/tlhsx/S-DCNN.

The validation methods in this paper were 5-fold cross-validation and independent testing. For the evaluation of the prediction results, we adopted the evaluation indicators commonly used in the identification of ATP binding residues: sensitivity (S _n ), specificity (S _p ), accuracy (ACC), and Matthews correlation coefficient (MCC) (i.e., Equations 8–11) (Zou et al., 2023). S n = T P T P + F N × 100 % (8) S P = T N T N + F P × 100 % (9) A CC = T P + T N T P + T N + F P + F N × 100 % (10) M C C = T P × T N − F P × F N T P + F P T P + F N T N + F P T N + F N (11) where, the number of ATP binding residues correctly predicted is TP, otherwise it is FN; the number of ATP non-binding residues correctly predicted is TN, otherwise it is FP. In addition, the flowchart was clearly described in Figure 4.

The basic feature parameters were input into the S-DCNN predictor, and the results of the 5-fold cross-validation were shown in Table 1. Here, the S_n and MCC values were 43.39% and 0.4101, respectively.

To improve the prediction performance, the dihedral angle, energy and propensity factors were introduced. The extracted dihedral angle, energy, propensity factor feature parameters and basic feature parameters were fused and input to the S-DCNN predictor, and the results were shown in Table 1 on 5-fold cross-validation.

In Table 1, when the feature parameters PP, E or F were respectively added to the feature B, the prediction results of B + F were relatively better. Then, the parameter PP or E was added to the above parameter set of B + F, and the results with parameter set of B + F + PP were relatively better. When the feature parameters PP, E and F were added at the same time, the best prediction results were obtained. The S_n, S_p, ACC and MCC values with feature set of B + F + PP + E reached 50.9%, 98.01%, 96.31% and 0.4773, respectively.

The prediction results of the three hyperparameters of the 5-fold cross-validation were shown in Figure 5. Figure 5A was a bar chart of the MCC and S_n values changing with the number of convolution layers. When the number of layers was 3, the S_n and MCC values reached the peak at the same time, then the optimal number of layers is 3. From Figures 5B, C, the optimal filters and batch size were 16 and 32, respectively.

The prediction results after optimization of hyperparameter were shown in (B + F + PP + E)* of Table 1. The S_n, S_p, ACC and MCC values reached 58.82%, 98.4%, 96.97%, and 0.5681, respectively.

To assess the efficacy of the SMOTE algorithm in predicting ATP binding residues, we conducted a comparative analysis between SMOTE and random undersampling alongside samples without preprocessing (referred to as RUS-DCNN and DCNN, respectively). In RUS-DCNN, for result stability, negative set samples were randomly selected ten times, with the average outcome of these ten selections serving as the final prediction. The prediction results of 5-fold cross-validation after optimization of hyperparameter were listed in Table 2. In Table 2, the MCC values of RUS-DCNN, DCNN and S-DCNN reached more than 0.438, and the ACC values of DCNN and S-DCNN reached more than 96.97%.

To verify the superiority of the S-DCNN algorithm, we computed the results of the SVM and random forest (RF) algorithm with SMOTE(i.e., S-SVM and S-RF) through 5-fold cross-validation, as detailed in Table 2. Specifically, the RF model utilized 500 decision trees, the SVM model employed a radial basis function kernel, and other parameters remained at default values.

To test the generalization ability of the prediction model, an independent testing set was utilized to predict ATP binding residues with the corresponding results outlined in Table 2 within brackets. Across the evaluation metrics of MCC and ACC, S-SVM, S-RF, and S-DCNN achieved values exceeding 0.409% and 96.2% respectively. Notably, S-DCNN demonstrated superior performance in terms of S_n and MCC. Moreover, the prediction model’s performance was assessed using the area under the Receiver Operating Characteristic (ROC) curve (AUC). Figure 6 illustrates the ROC curves for various algorithms based on SMOTE on the ATP-289 independent testing set, where S-SVM, S-RF, and S-DCNN yielded AUC values of 0.8585, 0.8841, and 0.9088 respectively.

To further verify the prediction performance of S-DCNN, we applied S-DCNN on another two frequently used datasets. The first set was constructed by Yu et al. (2013), in which the training set had 221 protein chains (ATP-221), and the independent testing set had 50 protein chains (ATP-50). The other set was constructed by Hu et al. (2018), in which the training set had 388 protein chains (ATP-388), and the independent testing set had 41 protein chains (ATP-41).

Using 5-fold cross-validation with optimized hyperparameters, the S-DCCN method was performed on the ATP-221 and ATP-388 datasets. The corresponding two prediction results were shown in Tables 3, 4, respectively. In Table 3, the S-DCNN achieved ACC of 97.0% on the ATP-221 dataset, surpassing other methods by 0.6%–0.8%. The MCC of the S-DCNN also exhibited significant improvement, with an increase ranging from 3.6% to 12.5%. In Table 4, the ACC and MCC values of the S-DCNN on the ATP-388 dataset reached 97.04% and 0.5887, respectively. To make a better comparison, we also listed the prediction results of the previous on the ATP-221 and ATP-388 datasets.

The prediction results of independent testing were listed in Tables 3, 4, with values displayed in brackets. Notably, Table 3 showcased the enhanced prediction performance of the S-DCNN method. Specifically, on the independent testing set ATP-50, S-DCNN achieved an ACC value of 97.0%, surpassing other methods by 0.2%–0.5%. Concurrently, the S_n and MCC values of S-DCNN exhibited notable enhancements, reaching 50.2% and 0.573, respectively, marking a 0.1%–6.5% increase and 3.9%–8.2% improvement compared to alternative methods. The S_p value of the S-DCNN was slightly higher than that of other methods. In Table 4, the ACC and MCC values of the S-DCNN on the independent testing set ATP-41 reached 96.78% and 0.5850 respectively. In addition, we drew the ROC curve of the S-DCNN method on the ATP-50 and ATP-41 sets, as shown in Figure 7. The AUC values of the S-DCNN method on the ATP-50 and ATP-41 sets were 0.9138 and 0.8973 respectively.