To validate the robustness and effectiveness of the proposed model, we performed it on three plant PPI datasets, A. thaliana, Z. mays, and O. sativa. The A. thaliana dataset was collected from TAIR¹ (Rhee et al., 2003), IntAct² (Kerrien et al., 2012), and BioGRID³ (Stark et al., 2006). After removing the redundancy, the final A. thaliana–positive dataset comprised 28,110 PPI pairs containing 7,437 A. thaliana proteins. These protein-interacted pairs constructed the primary A. thaliana PPI network. For the construction of the negative dataset, we employed a bipartite to formulate a network of plant PPIs, where the nodes represent the plant proteins and the links denote the interactions between them. Here, we use A. thaliana as an example. The whole associations between the 7,437 proteins are 55,308,969 (7,437 × 7,437) in the corresponding bipartite. However, only 28,110 PPIs had been demonstrated to have the interactions. Thus, the possible number of negative pairs is 55,280,859 (55,308,969–28,110), which is significantly more than the positive samples. To handle this binary classification problem, we randomly collected 28,110 non-interacting pairs as the negative dataset. In theoretical terms, the negative samples may contain a small number of positive samples; however, given the size of the whole non-interaction pairs, the probability of this situation is very small. In this way, the whole A. thaliana dataset consists of 56,220 protein pairs.

Maize and rice are the main cash crops in the world. The maize (Z. mays) dataset contains 14,800 positive pairs, which was downloaded from PPIM⁴ (Zhu et al., 2016) and agriGO⁵ (Tian et al., 2017). Similarly, we assumed that the proteins in different subcellular work compartments have no interactions and finally achieved 14,800 non-interacting protein pairs. The rice (O. sativa) dataset consisted of 9,600 protein pairs, 4,800 positive pairs, and 4,800 negative pairs collected from the PRIN database⁶ (Gu et al., 2011).

To mine highly efficient features for training the models, each protein pair is encoded as 800-dimensional feature vector by PSSM (Gribskov et al., 1987). PSSM has been successfully employed in various fields of biological research including the prediction of PPI site, subcellular localization, and DNA-binding protein identification. In this section, we applied PSI-BLAST (Altschul and Koonin, 1998) tool to represent protein sequence as a U × 20 matrix, where Q = {η_a,b:a = 1⋯U andb = 1⋯20}, and it can obtain the information of plant sequential evolution. PSSM can be defined as

where η_a,b represents probability that the a-th mutate to b-th amino acid during the evolutionary process. In the experiment, plant protein sequences were adopted as seeds to search and align homogenous sequences from SwissProt database by PSI-BLAST tool. The tool will be used to recognize members of gene family and evolutionary relationships between plant protein sequences. It is also able to generate a 20-dimensional vector to denote the probabilities of conservation against mutations to the 20 amino acids. The number of iterations is set to 3 and the E-value is cut off at 0.001 to achieve homologous sequences. The PSI-BLAST tool and SwissProt database can be accessed online⁷.

In this section, we introduce discrete Hilbert transform (DHT; Cizek, 1970) to extract feature descriptors from the PSSM to make the prediction more convenient and accurate. DHT is used as a tool for signal analysis in the time and frequency domains. Before describing the 2-dimensional DHT, the 1-D DHT (Ponomareva et al., 2018) is used in the spatial and frequency domain and has been previously described (Stark, 1971; Bracewell and Bracewell, 1986; Zhu et al., 1990; Onodera et al., 2005).

To better extract the feature descriptors, we used the 2-D DHT for constructing the local energy of PSSM. In this work, we applied the 2-D DHT, which is defined by Read and Treitel (1973) in the frequency domain. Our Matlab code is shown as follows:

function x = hilbert2(xr,m,n)

%HILBERT2 Discrete-time 2D analytic signal via Hilbert transform.

% X = HILBERT2(Xr) computes the 2D discrete-time analytic signal

% X = Xr + i^∗Xi such that Xi is the Hilbert transform of real image Xr.

% If the input Xr is complex, then only the real part is used: Xr = real(Xr).

% HILBERT2(Xr,M,N) computes the MxN-point Hilbert transform. Xr is padded

% zeros if it has less than MxN points, and truncated if it has more.

if nargin < 2, n = []; end

if ∼isreal (xr)

   warning (’HILBERT2 ignores imaginary part of input.’)

   xr = real (xr);

end

if isempty (n)

   [m, n] = size (xr);

end

if m < 2 | | n < 2,

   x = Hilbert (xr); % 1D analytic signal

   return;

end;

In this work, PSI-BLAST encoded each protein sequence as a U×20 matrix. Due to the different lengths of protein sequences, the size of each matrix constructed by PSSM is also different. To handle this problem, we transformed the variably sized PSSM into a 20×20 matrix, and the 2-D DHT is applied to extract feature vectors from the PSSM profile. In this way, each plant protein sequence will be converted into a 400-dimensional vector by 2-D DHT. As a non-linear filtering technique, SVD has been widely applied in noise reduction of vibration signals. This is because the signals after noise reduction have a small phase-shift and there is no time delay effect. To improve the prediction accuracy and reduce the dimensionality of the input feature matrix, we applied SVD (Klema and Laub, 1980) algorithm to reduce size of feature vectors from 400 to 300. At the same time, the lower dimensions could reduce the complexity of the model and increase the generalization error of the classifier. Finally, each protein pair will be represented as a 600-dimensional DHT descriptor.

Considering the larger numbers of hidden layers that can be used for training networks, artificial neural networks consist of two or more hidden layers that are often referred as DNN as shown in Figure 1. The depth of a neural network relates to the quantity of hidden layers, and the largest number of neurons determines the width of DNN (Hinton et al., 2006; Hinton and Salakhutdinov, 2006).

In terms of structure, DNN is composed of many plain modules, which appear as a multilayer stack. The data are first received by the input layer, and then converted through a non-linear way across many hidden layers. Before calculating the final output, the average gradient is first computed and the corresponding weights are adjusted. Neurons of a hidden layer or input layer are associated with the neurons of the existing layer. Each neuron will compute a weighted sum of its input and perform a non-linear activation function to capture its outputs. The non-linear activation functions usually include sigmoid, rectified linear unit (ReLU), and hyperbolic tangent. In this work, we used the sigmoid and ReLU. We constructed a DNN-based model using the TensorFlow platform shown in Figure 1. This model consists of two hidden layers with 48 neurons each. The DHT feature descriptors are employed as the inputs for the DNN model. After that, these features were set into the hidden layers for training and predicting PPIs. Adam algorithm (Kingma and Ba, 2014), which is an adaptive learning rate approach, was adopted in our methods to accelerate the training process. At the same time, to avoid overfitting, the dropout technique was also applied to our model (Khan et al., 2019). We also used the cross-entropy loss and ReLU activation function to speed our training and achieve better predictive performance (Hinton et al., 2015). The loss can be calculated by the following formulas:

In Eqs. 2–6, n describes the amount of protein pairs that need to be trained, m denotes the individual network, h₁ represents the depth of two individual networks, and h₂ denotes the depth of the fused network. The activation function of ReLU and output layer with sigmoid is σ₁ and σ₂, respectively; ⊕ is the concatenation operator. R represents the output of hidden layer and y is the corresponding desired output. T and b indicate the weight matrix and bias vectors.

To prevent overfitting and validate the robustness of our method, five-fold cross-validation (CV) scheme is performed on our method. Specifically, the entire plant’s PPI dataset will be randomly split into five equal parts; four of them will be employed for training and the remaining one was used for testing. The training and testing data will not overlap with each other to prevent overfitting. The final validation results were the mean value obtained by the five-fold CV scheme. The predictive performance of the proposed approach is verified by five different measurements, including accuracy (Acc), precision (PR), sensitivity (Sens), specificity (Spec), and MCC. They can be represented by

where TP, FP, TN, and FN are associated with the number of true positive, false negative, true negative, and false negative, respectively. In addition, receiver operating characteristic (ROC) curves (Hand, 2009) were plotted for better accessing the predictive performance of the proposed model. Furthermore, AUC (area under ROC curve) Huang and Ling (2005) values were also used as an evaluation criterion.

We validated the predictive performance of the proposed model on three plant PPI datasets by five-fold CV scheme, including A. thaliana, Z. mays, and O. sativa. It can be observed from Table 1 that the average accuracy (Acc), precision (PR), sensitivity (Sens), specificity (Spec), and Matthews correlation coefficient (MCC) and AUC values obtained on the A. thaliana dataset are 71.48%, 66.64%, 86.09%, 56.88%, 44.94%, and 0.8369, respectively. Their SDs are 0.69, 0.89, 1.08, 2.21, 1.14, and 0.36%, respectively. Table 2 lists the prediction results obtained on the Z. mays dataset, from which we can see the average Acc of 85.41%, PR of 81.54%, Sens of 91.67%, Spec of 79.17%, MCC of 71.43%, and AUC of 0.9466, respectively. Their SDs are 1.18, 2.38, 1.18, 3.24, 2.00, and 0.26%, respectively. On the O. sativa dataset, shown in Table 3, our model performs at an Acc of 82.60%, PR of 75.79%, Sens of 95.89%, Spec of 69.31%, MCC of 67.65%, and AUC of 0.9442, with SDs of 1.79, 2.43, 0.91, 3.53, 2.98, and 0.58%, respectively. Figures 2–4 illustrate the ROC curves yielded on A. thaliana, Z. mays, and O. sativa datasets. In the figure of ROC curves, x-axis is the false positive rate and y-axis represents the true positive rate.

Based on the experimental results, it can be indicated that the proposed model is effective for identifying PPIs in plants. We attributed this better prediction performance to the powerful DHT–SVD descriptors and the excellent DNN classifier. The PSSM not only encodes the sequence into matrix but also obtains the sufficient prior information of plant proteins. In addition, the application of DHT extracted robust feature descriptors from PSSM, and then, SVD algorithm was employed to reduce the noise and decrease the dimension of feature matrix that can better improve the prediction performance. As a popular deep learning classifier, DNN shows the powerful ability for training and predicting, which makes us more convinced that our method can be a useful tool for plant PPI prediction.

There are many machine learning classifiers that have been applied to predict PPIs. K-nearest neighbor (KNN) (Keller et al., 1985) and random forest (RF) (Breiman, 2001) are the most widely used algorithms. The KNN algorithm is one of the simplest classification approaches and it has been widely applied to detect PPIs (Li et al., 2009). RF is a decision tree–based ensemble learning method, and it is known for its powerful ability of classification (Li et al., 2012). To further verify the predictive ability of DNN classifier, we compared it with the KNN and RF model by the five-fold CV scheme and adopted the same DHT feature descriptors. The results list in Table 4 illustrates that our method achieved higher AUC values across the A. thaliana, Z. mays, and O. sativa datasets. It can be observed that the average AUC values of the DNN classifier are 0.1023, 0.1215, and 0.1354 higher than those of KNN classifier. Similarly, when compared with the RF classifier, the AUC value of our model improved 0.0036, 0.013, and 0.0241, respectively. From the comparison results shown in Figure 5, we considered that the combination of DNN classifier and DHT descriptors can significantly improve the performance in plant PPI prediction.

To evaluate the performance of PSSM, we compared it with the substitution matrix representation (SMR), which was proposed by Yu et al. (2012) to represent protein sequence. In this section, we employed the BLOSUM62 matrix to encode the A. thaliana protein sequence as a 20 × 20 matrix. Then, the DHT algorithm was applied to extract feature descriptors from SMR matrix and SVD was also adopted to reduce the feature dimensions. By this way, we can generate a 600-dimensional SMR–DHT descriptor for each protein pair. The five-fold CV results of SMR–DHT descriptors combined with DNN classifier on the A. thaliana dataset are summarized in Table 5. It can be observed that the PSSM-based method performs significantly better than the SMR-based method. For example, the accuracy and AUC gaps between PSSM and SMR-based method are 4.38 and 4.94%, respectively. The higher predictive accuracy and lower SDs further indicated that our method performs better than the SMR-based approach (Figure 6).

To illustrate the effectivity of our feature extraction approach, we compared DHT with some popular correlative methods, including discrete cosine transform (DCT) (Ahmed et al., 1974), fast Fourier transform (FFT) (Nussbaumer, 1981), discrete wavelet transform (DWT) (Nanni et al., 2012), and auto-covariance (AC) (Zeng et al., 2009). As shown in Table 6 and Figure 7, on the O. sativa dataset, our method obtained a high prediction accuracy of 82.60%. The prediction accuracy values of other methods are 80.95, 75.31, 81.54, and 66.63%, respectively. Our method performs better than the other four methods. Especially compared with the AC-based method, our approach improved the Acc, Spec, MCC, and AUC by over 15%, and PR and Sens by over 7%, respectively. Although the Sens value of our method is not the highest, it still obtains an excellent value of 95.89%. The Acc, PR, Sens, Spec, MCC, and AUC values obtained from our model are 1.06, 0.69, 1.08, 1.05, 2.15, and 1.31% higher than the values of the DWT-based method. These comparison results further indicated the superiority of the proposed method.

To further validate the potential of the presented method, we performed it on the yeast and human PPI dataset, which was introduced by Guo et al. (2008) and Huang et al. (2015). The predictive results of the two datasets are listed in Tables 7, 8, and the corresponding ROC curves are shown in Figures 8, 9. When performing on the yeast dataset, it achieved average Acc, PR, Sens, Spec, MCC, and AUC value of 79.54%, 73.46%, 92.63%, 66.47%, 61.27%, and 0.9203, with SDs of 1.43, 2.11, 1.18, 3.60, 2.16, and 0.46%, respectively. From Table 8, it can be observed that the proposed model yielded great results on the human dataset, an average Acc of 82.76%, PR of 75.79%, Sens of 94.18%, Spec of 72.30%, MCC of 67.72%, and AUC of 0.9473, with SDs of 1.68, 2.79, 1.64, 4.74, 2.37, and 0.26%, respectively. From these results, we can observe that the powerful DNN-based classifier combined with the DHT feature descriptor is accurate and robust for exploring cross-species predictions of PPIs.