Computer-aided prediction of antigen presenting cell modulators for designing peptide-based vaccine adjuvants.

BACKGROUND: Evidences in literature strongly advocate the potential of immunomodulatory peptides for use as vaccine adjuvants. All the mechanisms of vaccine adjuvants ensuing immunostimulatory effects directly or indirectly stimulate antigen presenting cells (APCs). While numerous methods have been developed in the past for predicting B cell and T-cell epitopes; no method is available for predicting the peptides that can modulate the APCs.
METHODS: We named the peptides that can activate APCs as A-cell epitopes and developed methods for their prediction in this study. A dataset of experimentally validated A-cell epitopes was collected and compiled from various resources. To predict A-cell epitopes, we developed support vector machine-based machine learning models using different sequence-based features.
RESULTS: A hybrid model developed on a combination of sequence-based features (dipeptide composition and motif occurrence), achieved the highest accuracy of 95.71% with Matthews correlation coefficient (MCC) value of 0.91 on the training dataset. We also evaluated the hybrid models on an independent dataset and achieved a comparable accuracy of 95.00% with MCC 0.90.
CONCLUSION: The models developed in this study were implemented in a web-based platform VaxinPAD to predict and design immunomodulatory peptides or A-cell epitopes. This web server available at http://webs.iiitd.edu.in/raghava/vaxinpad/ will facilitate researchers in designing peptide-based vaccine adjuvants.

Background

Peptide subunit vaccines are hailed as an advancement over live or inactivated whole organism vaccines due to their ability to minimize adverse reactions [1]. Yet, antigenic peptides by themselves are poorly immunogenic since they lack the capability of activating the innate immunity. Activation of the innate immune system is required for stimulation of whole immune system including adaptive immunity. Hence, there is a need for inclusion of immunostimulants known as adjuvants in the subunit vaccine formulations. Conventionally, empirical approaches were used for adjuvant discovery, so far limited adjuvants have been approved and licensed for clinical use like alum, MF59, AS03 and AS04 [2].

Vaccine adjuvants effectuate their action by a variety of mechanisms with all of them involving the antigen presenting cells (APCs) particularly the dendritic cells [3]. One of these mechanisms is the activation of the pattern recognition receptors (PRRs) on the APCs that recognize conserved microbial molecular signatures. PRR ligands shape the adaptive immune response mediated by the APCs. A majority of the vaccine adjuvants are ligands of PRRs making them potential targets for rational design of vaccine adjuvants [2]. Thus, hypothesis-driven adjuvant development relies on the expectation that the mechanistic understanding of the immune responses exhibited by PRR ligands would enable fine-tuning the specificity of adjuvants to attain vaccine efficacy and safety, simultaneously. An important example of a class of molecules that have been shown to have immunomodulatory effects and are poised to become safe and cost-effective adjuvants in future is—short immunomodulating peptides [4]. Figure 1 is a schematic representation of the adaptive immune cell activation by a coordination of antigen presentation to the naïve adaptive immune cell with the release of cytokine milieu mediated by PRR activation. Keeping in view the role of peptide ligands of PRRs in the activation of APCs, we introduce the term ‘A-cell epitopes’ for these immunomodulatory peptides.

Fig. 1

An illustrative mechanism of antigen presenting cell (APC) activation caused by immunomodulatory peptides through innate immune receptors leading to the induction of adaptive immune cells. The immunomodulatory peptides are ligands of innate immune receptors that evoke cytokine expression through cellular signaling pathways. The cytokines lead to the maturation of naïve cells into mature adaptive immune cells such as various types of T-lymphocytes. Since the immunomodulatory peptides activate the APCs leading to the activation of the adaptive immune cells, they may be used as vaccine adjuvants and be called ‘A-cell epitopes’. The figure was drawn using ScienceSlides, made available at http://www.scienceslides.com/ by VisiScience

Cationic host defense peptides (HDPs) were originally discovered as antimicrobial peptides produced within the multicellular organisms having a broad-spectrum activity against bacteria, viruses, fungi, protozoa, etc. [5]. Of late, HDPs and their synthetic analogs called innate defense regulators (IDRs) have been realized to cause immunomodulatory effects like differentiation and activation of innate and adaptive immune cells, modulation of pro- and anti-inflammatory responses, chemo attraction, autophagy, apoptosis and enhancement of immune-mediated bacterial killing [6]. Many host defense peptides (HDPs) with known immunomodulatory effects are already in clinical trials [7]. IMMACCEL-R is a short synthetic peptide with immunomodulatory properties that has been commercialized for use as vaccine adjuvant in animals and birds for the purpose of antibody generation [8]. The human cathelicidin antimicrobial peptide (CAMP or LL37) is another well-known antimicrobial peptide shown to induce immunomodulatory effects [9] and has been found to be associated with immune-related disorders like psoriasis [10] and morbus Kostmann [11].

Adjuvants have been incorporated into the vaccine formulations for qualitative alteration of the adaptive immune responses that are different from non-adjuvanted antigens. For instance, the adjuvants have been used to skew the immune responses with respect to Th1 (T-helper 1) cells versus Th2 cells, CD8+ versus CD4+ T cells, specific antibody types, etc. [2]. The PRR ligands produce these effects by virtue of the adapter proteins in the signaling pathways activated by the PRR. For example, bacterial flagellin protein causes adjuvant effect through TLR 5 and produces a mixed Th1 and Th2 response instead of polarized Th1 response and requires TLR adapter protein MyD88 for this effect. In contrast, monophosphoryl lipid A (MPL) and bacterial lipopolysaccharide (LPS) acting through TLR 4 activation lead to production of pro-inflammatory cytokine TNF leading to a polarized Th1 response instead of mixed Th1–Th2 response. While MPL signals through TRIF adaptor, LPS mediated activation of TLR 4 acts through both TRIF and MyD88 adapter proteins. Thus, MPL formulated on alum (AS04) stimulates a polarized Th1 cell response and is a component of licensed vaccine for HBV and papilloma that has proven to be both safe and effective.

Developing adjuvants based merely on empirical studies without the understanding of mechanisms is inadequate [2]. There is a need to develop systematic and rational approaches for designing highly potent vaccine adjuvants. One such approach could be the development of PRR ligands into vaccine adjuvants since their mechanism is known.

In such a scenario, in silico models to screen and identify potential vaccine adjuvant candidates could prove to be useful as the existing experimental approaches are time and resources consuming [12]. Previously, our group developed a method, Vaccine DA for predicting immunomodulatory oligodeoxynucleotides that can activate innate immune system via Toll-like receptor-9 (TLR-9). This tool can be used for designing oligodeoxynucleotide-based vaccine adjuvants as well as for genome-wide screening of vaccine adjuvants [13]. Recently, we also developed a method imRNA for designing single-stranded RNA (ssRNA) based vaccine adjuvants [14]. These methods may play an important role in designing DNA and RNA-based therapeutics as these methods allow a user to design oligonucleotides of desired immunogenicity.

In the last two decades, numerous method have been developed for predicting potential of peptides to stimulate adaptive arm of the immune system that include methods for predicting MHC binders, B-cell epitopes [15-23] and T-cell epitopes [24-31]. To the best of our knowledge, no method has been developed so far for predicting immunostimulatory potential of peptides to activate innate immunity. In this study, we made an effort to develop method for predicting immunomodulatory peptides that can activate innate arm of immune system or antigen presenting cells. These peptides activate the antigen presenting cells (e.g., dendritic, macrophages); hence, we propose that these immunomodulatory peptides be termed as ‘A-cell epitopes’.

In the present work, first we collected experimentally identified immunomodulatory peptides from the literature and included them in our positive set named A-cell epitopes. Next, we collected the human endogenously circulating peptides to build the negative set named A-cell non-epitopes. Combining the positive and the negative sets into a complete dataset, we developed support vector machine (SVM) based computational models that can classify a new query peptide as A-cell epitope or non-epitope. To benefit the users of the scientific community, we provided the best performing SVM-based prediction models in the form of a web-based application called VaxinPAD to be used for identifying and designing novel A-cell epitopes. Such peptides identified computationally might serve as the starting molecules for designing peptide-based vaccine adjuvants.

Methods

Dataset

The experimentally validated immunomodulatory peptide sequences were obtained from 16 patents. As an example of the sequences considered immunomodulatory, a set of sequences taken from a patent (US20110008318 A1), includes flagellin-derived peptides that exhibit immunomodulatory effect by direct binding to TLR 5 as indicated by assays reporting increased NF-κB expression estimated from coupled luciferase activity and TNFα production estimated using flow cytometry. In one of the patents (US7462360 B2), a class of immunomodulatory peptides, called alloferons, derived from the bacteria challenged blood of larvae of the insect blowfly, Calliphora vicina, have been found to stimulate the cytotoxic anticancer activity of the human NK-cells and lymphocytes. In another case, a set of peptides as described in patent US8791061 B2 have been shown to enhance innate immunity by modulating the activity of type II transmembrane serine protease dipeptidyl peptidase (DPPIV) also known as CD26 or adenosine deaminase binding protein, expressed on major immune cells like activated T-cells, B-cells, NK-cells, macrophages and epithelial cells. With two major functions of signal transduction and proteolysis, the effects of DPPIV protein-mediated cellular processes include modulation of the chemokine activity by cleaving dipeptides from chemokine N-terminus that alters the receptor binding and specificity of the processed chemokine. DPPIV is a neutrophil chemorepellant and eosinophil chemoattractant too.

After removing the longer sequences, 304 unique sequences left in the length range of 3–30 residues were used to constitute the positive dataset named here as the A-cell epitopes. The upper bound of length 30 residues was kept as more than 90% of the originally collected epitope sequences were retained keeping this criterion used for removing very long sequences. In the absence of experimentally verified non-immunomodulatory peptides (non-epitopes), the experimentally identified endogenous human serum peptides [32, 33] were taken as non-epitopes. We assume these peptides are non-immunogenic as they are part of human serum, thus we assign them as non-epitopes. Only the sequences of the length 3–30 were taken into the negative dataset. In this manner, the main dataset consisted of 304 A-cell epitopes and 385 non-epitopes. Additional file 1: Table S1 provides the sequences and the source patent/publication for the positive and the negative datasets.

Input features

In order to develop any in silico model it is important to generate input features corresponding to each data point. In this study, a data point is the amino acid sequence of a peptide (either A-cell epitope or non-epitope). It is important to generate fixed length input features because machine-learning techniques require fixed length vector for developing a model. As the length of peptides is variable, thus we computed amino acid composition of A-cell epitopes and non-epitopes for developing models. We also computed the average amino acid composition of A-cell epitopes, non-epitopes and the human proteins, in order to understand compositional bias in A-cell epitopes. The amino acid composition for each sequence constituted the input vector of length 20, which was used for developing SVM-based prediction models. Similarly, the dipeptide composition vectors of length 400 were generated for A-cell epitopes and non-epitopes with each element of a vector corresponding to the composition value of each type of possible dipeptide. In addition to compositional features, we also generated binary features for developing models using fixed length of amino acids from the termini (N-terminal or C-terminal or both) of peptides. In the case of binary feature, an amino acid is represented by a vector of 20, where the presence of amino acid is indicated by ‘1’ and the absence is presented by ‘0’ [34]. This means a peptide of length N is presented by a vector of length N × 20 in the case of binary features.

Motif search

We used the Motif—EmeRging and with Classes—Identification (MERCI) Program [35] to identify motifs exclusively occurring in the A-cell epitopes [36]. Though this program allows searching for gapped and ungapped motifs, but we restricted our analysis to the ungapped motifs. It is well established that in the case of T-cell epitopes, even a single residue mutation changes its immunogenicity [37] and can even eliminate the immunogenicity of the epitope [38]. Hence, intuitively the ungapped motifs found to be conserved among the positive sequences are more likely to help identify novel A-cell epitopes. Thus, we computed and compared the frequency of occurrence of ungapped motifs in A-cell epitopes, non-epitopes and the Swiss-Prot proteins.

Classifiers based on machine learning techniques

In the present study, some commonly used popular machine learning techniques were used to develop classification prediction models. We used WEKA package to implement these machine learning techniques namely Random Forest, Naïve Bayes, SMO and J48 [39]. These classification models were developed using commonly used features of peptides like amino acid composition (AAC) and the dipeptide composition (DPC).

Support vector machine (SVM)

Subsequent prediction models in this study were developed using SVM, which has been frequently used to develop models for epitope prediction in previous studies [21, 23, 28]. SVM has been the method of choice for building epitope prediction models especially T-cell epitopes [40] due to its ability to provide effective models on high dimensionality data with less data points. Also, in the past studies it has been shown that SVM performs better on independent dataset in comparison to other machine learning classifiers [41]. The dataset used in the current study contains data points comparable in number to the dimensionality. Hence, we optimized the prediction models on various parameters using the radial basis kernel of a freely available program SVM [42] to select the best performing models on different sets of features.

Evaluation of models using internal and external validation

In this study, standard procedure was followed to evaluate the performance of models in order to avoid biases in performance due to over optimization. Our main dataset was divided into two categories internal and external dataset, where the internal dataset contained ~ 80% sequences and the external dataset comprised of the remaining 20% sequences. In order to perform internal validation, we performed fivefold cross validation technique on internal dataset. In this technique, the dataset is divided in five sets, four sets are used for training a model, and the remaining set is used for testing the model. This process is repeated five times so each sequence is tested only one time. In order to perform the external validation of a model, the best model developed using fivefold cross validation is tested on an external dataset. It is important to assess the performance of a model on external or independent dataset because the performance of a model in internal validation may be biased due to optimization of the model [28]. The performance of models was measured using standard threshold dependent parameters namely sensitivity, specificity, accuracy and Matthew’s correlation coefficient (MCC) [19, 36] and a threshold independent parameter area under receiver operating characteristics (AUROC) [43].

Bootstrap aggregating

In order to avoid over fitting of models and reducing variance in performance of models; we used bootstrap aggregating (bagging) for averaging performance of models. In this study, process of creating internal and the external datasets has been repeated ten times. Each time, the sequences for the internal dataset were randomly selected from the main dataset, and the remaining sequences were included in external dataset. Finally, we evaluated the performance of our models using various features on both the internal as well as the external datasets as described in above sections. This process gave 10 performance values using internal and 10 performance values using external validation from 10 rounds of sampling. We computed the mean and standard deviations of these performance values to check for bias in performance of the models depending on the choice of sequences on which the models were trained or independently evaluated.

Random peptides as negative dataset

As described above, initially the negative dataset consisted of the experimentally identified endogenous human serum peptides as non-epitopes constituting the negative dataset. We further wanted to check whether the performances of the classification models were dependent on the choice and size of the negative datasets. This was necessary as the negative dataset does not contain the experimentally verified non-epitopes. For this, we created an alternative negative dataset of random peptides derived from the human proteins obtained from the Swiss-Prot database. As mentioned in the previous section, for each of the 10 rounds of sampling, a different set of random peptides 10 times the number of the positive sequences (A-cell epitopes) from the human proteins was kept as the negative dataset.

Results

Compositional analysis

One of the objectives of this study is to understand the nature of A-cell epitopes regarding the residues preferred in A-cell epitopes. Thus, we computed the average residue composition of A-cell epitopes and the non-epitopes. The non-epitope dataset consists of peptides occurring in the normal human serum assumed to be non-immunomodulatory. In addition, the average residue composition of the Swiss-Prot Human proteins was also computed and compared with that of the A-cell epitopes.

In the A-cell epitope dataset, the percentage composition of an amino acid residue was calculated for each epitope, and the average of these values was plotted in Fig. 2 for the corresponding amino acid. Similarly, the average percentage composition was calculated for all the amino acids in the non-epitope dataset and the Swiss-Prot Human proteins. As shown in Fig. 2, the residues showing noticeable differences in average composition between A-cell epitopes and non-epitopes are C, D, E, I, L, R, S, T, V and W. Student’s t-test significance value (p-value) was calculated for each residue type to check whether the composition values among A-cell epitopes were different from those in the non-epitopes. In decreasing order of significance (increasing adjusted p-value), the residues R, E, T, S, D, V, W, L, I and C showed the most significant difference between the A-cell epitopes and non-epitopes among all of the residue types with adjusted p-values 1.76E−39, 1.66E−24, 7.04E−18, 3.59E−17, 3.72E−12, 2.79E−11, 5.21E−10, 7.91E−10, 1.24E−09, 1.76E−08 respectively (Additional file 1: Table S2). In particular, when compared to the human proteins taken from Swiss-Prot; R was found to have a higher average composition in A-cell epitopes. The average composition of R in non-epitopes is lower than that in the human proteins. Overall, the residues I, R, V and W were found to be more abundant in the A-cell epitopes as compared to the non-epitopes and Swiss-Prot Human proteins.

Fig. 2

Barplots showing the comparison of percent average amino acid composition of A-cell epitopes (blue) with non-epitopes (red) and Swiss-Prot human proteins (green)

Similarly, the dipeptide and tripeptide compositions of the A-cell epitopes and non-epitopes were also compared with the Swiss-Prot Human proteins. Additional file 1: Table S3 gives the average composition for each dipeptide in the A-cell epitopes, non-epitopes as well as the Swiss-Prot Human proteins. After sorting the table according to descending order of difference of dipeptide composition between the A-cell epitopes and the human proteins, top 10 dinucleotide include the residues I, R and V. But these motifs also contain other amino acids that show less significant difference of abundance as compared to the non-epitopes and human proteins. Similar analysis of tripeptide composition is shown in Additional file 1: Table S4. In this case too, the top 10 tripeptide motifs include less abundant residues apart from I, R and V.

Terminal residue preference

We performed position-specific analysis of residues in A-cell epitopes to understand the type of residues preferred at different positions in A-cell epitopes. In this study, two-sample logo (TSL) tool (available at http://www.twosamplelogo.org/cgi-bin/tsl/tsl.cgi) [44] was used to visualize residues preferred or not preferred in A-cell epitopes. Since the minimum peptide length in the dataset was 3, the N-terminal 3 residues of both the negative and the positive sequences were taken as input to build the N-terminus TSL. C-terminus TSL was obtained using the C-terminal 3 residues from the dataset. Figure 3 shows that the residues R, V and I are among the preferred residues in the A-cell epitopes at both the N and the C termini.

Fig. 3

Two-sample logo of the 3 residue positions at the a N-terminus and b C-terminus of the A-cell epitopes and non-epitopes. Enriched label represents the positive dataset whereas depleted label represents the negative dataset. In a two-sample logo, the height of a symbol at a residue position is proportional to the difference in symbol frequency between the positive and the negative datasets at that residue position. In the case of A-cell epitopes (as positives) and non-epitopes (as negatives), R is a preferred amino acid at terminal positions apart from I and V. The symbol colors are in accordance with the WebLogo default color scheme provided by the server available at http://www.twosamplelogo.org/cgi-bin/tsl/tsl.cgi. In the default WebLogo color scheme, residues G, S, T, Y and C appear in green color, N and Q are colored purple, K, R and H are depicted in blue, D and E are drawn red and P, A, W, F, L, I, M and V are shaded black

MERCI motif analysis

The Motif—EmeRging and with Classes—Identification (MERCI) Program is a software that helps in finding the motifs exclusive to one class when compared to another class of sequences. Additional file 1: Table S5 provides the MERCI motifs exclusive to the A-cell epitopes as compared to non-epitopes. Top 10 ungapped motifs with respect to the occurrence in the A-cell epitope sequences have a frequent occurrence of I, R and V. On the other hand, ungapped MERCI motifs exclusive in non-epitopes (Additional file 1: Table S6) that are top 10 in abundance contain E, G, P and L.

Rare motif occurrence

We compared the occurrence of peptide n-mers (n = 3, 4, 5, 6) in the A-cell epitopes and non-epitopes. First, the occurrence of each type of n-mer was counted in all of the Swiss-Prot proteins, and the n-mers were arranged in increasing order of occurrence. In this order, the n-mers were divided into 8 bins such that the 1st bin contained the n-mers least abundant in Swiss-Prot while the 8th bin contained the most abundant n-mers occurring in the Swiss-Prot. Next, the percentage of n-mers in a particular bin that occur in the Swiss-Prot was calculated with respect to the total number of n-mers in Swiss-Prot. Similar percentage value was calculated for A-cell epitopes and non-epitopes for each bin, and the values were presented in the form of a plot in Fig. 4.

Fig. 4

Comparison of occurrence of a tripeptides, b tetrapeptides, c pentapeptides and d hexapeptides divided into 8 bins in the ascending order of occurrence (most rarely occurring to most abundant) in Swiss-Prot proteins

Figure 4a shows that tripeptides in the first three bins occur more in A-cell epitopes while those of 4th, 5th, 6th, 7th and 8th bin (most abundant Swiss-Prot tripeptides) occur more in the non-epitopes. For tetrapeptides (Fig. 4b), the bins having more number of tetrapeptides occurring in the A-cell epitopes than non-epitopes are 1st, 2nd, 3rd and 4th. Figure 4c shows the occurrence of the pentapeptides. The bins having distinctly more pentapeptides in the A-cell epitopes than non-epitopes are again the first four bins. On the other hand, the percentage occurrence of hexapeptides of A-cell epitopes is lower than non-epitopes and Swiss-Prot proteins only in the 8th bin (Fig. 4d).

Prediction of immunomodulatory peptides

The sequence-based analyses like residue composition preferences; position-wise residue preference and motif search indicated that these features could help in discriminating the A-cell epitopes from non-epitopes. We developed SVM-based prediction models using SVM by from the dataset of 304 A-cell epitopes as positive sequences and 385 non-epitopes as negative sequences. From each of the positive and negative datasets, ~ 80% sequences were kept in the training–testing dataset while the remaining ~ 20% were kept in the independent dataset. Thus, the training–testing dataset had 243 positive and 308 negative sequences. The best performing models were selected on the basis of highest Matthews correlation coefficient values and a minimal difference between the sensitivity and specificity values.

Prediction models based on machine learning techniques

In order to understand, which machine learning technique will be most efficient for predicting A-cell epitopes, models were developed using different machine learning techniques. Initially models were developed using SVM implemented with SVMlight and four commonly used techniques (Random Forest, Naïve Bayes, SMO and J48) implemented using WEKA package. These models were developed using amino acid composition (AAC) and dipeptide composition (DPC) of peptide sequences (epitope and non-epitope). As shown in Additional file 1: Table S7, SVM based model performed better than models developed using any other machine learning technique. SVM based models on training dataset obtained MCC values 0.90 and 0.91 for AAC and DPC respectively. Similarly, performance was evaluated on the independent dataset. Thus, in this study, we used SVM for developing models using various features of peptides. The performance of SVM models developed using different features have been shown in Additional file 1: Table S8. A variation was observed in performances of models on some feature sets between the training and the independent datasets. Ideally models performing better on the training dataset should also perform better on independent dataset. In our case, we observed that models performing better on training dataset were performing lower on independent dataset and vice versa (see Additional file 1: Table S8). This inconsistency in models’ performance might have arisen from the over fitting of models to the training dataset.

Diminishing over fitting using bagging

In this study, we used bagging approach for sampling to overcome problem of over optimization or over fitting. Bagging procedure of sampling was adopted in this study to evaluate performance of models. The main dataset was divided ten times randomly into the internal and external datasets. Thus, we get ten training/internal datasets and ten independent datasets. The performance of all models was evaluated ten times on training and independent dataset. Finally, we computed the performance of each model in 10 rounds of sampling and reported the mean and standard deviation on these 10 datasets. Table 1 provides the performance values of models developed on various feature sets for the 10 internal dataset samples in different categories.

Table 1

The performance of SVM-based models developed using various features; models were evaluated on training dataset using fivefold cross-validation (internal cross-validation)

Features	Threshold	Sensitivity (%)	Specificity (%)	Accuracy (%)	MCC	AUROC	Parameters
AAC	− 0.1	94.49 ± 0.80	92.38 ± 1.33	93.30 ± 0.84	0.87 ± 0.01	0.98 ± 0.00	g = 0.001, c = 3, j = 3
N5 AAC	0	88.54 ± 0.75	90.25 ± 1.87	89.44 ± 1.26	0.79 ± 0.02	0.94 ± 0.00	g = 0.0005, c = 2, j = 1
C5 AAC	0	91.13 ± 1.42	92.94 ± 1.20	92.08 ± 1.05	0.84 ± 0.02	0.97 ± 0.00	g = 0.001, c = 9, j = 1
N5C5 AAC	− 0.2	93.73 ± 0.60	92.83 ± 0.76	93.26 ± 0.40	0.87 ± 0.00	0.98 ± 0.00	g = 0.0005, c = 1, j = 1
DPC	0	93.79 ± 1.12	95.68 ± 0.78	94.84 ± 0.72	0.90 ± 0.01	0.99 ± 0.01	g = 0.0005, c = 1, j = 2
N5 DPC	− 0.1	83.42 ± 1.77	87.73 ± 2.00	85.69 ± 1.10	0.71 ± 0.02	0.93 ± 0.00	g = 1e−05, c = 9, j = 1
C5 DPC	− 0.1	90.21 ± 0.91	93.62 ± 0.96	92.00 ± 0.50	0.84 ± 0.01	0.97 ± 0.00	g = 0.0005, c = 1, j = 2
N5C5 DPC	− 0.2	93.60 ± 0.72	92.67 ± 1.16	93.11 ± 0.70	0.86 ± 0.01	0.98 ± 0.00	g = 0.0001, c = 1, j = 1
N5 bin	− 0.1	86.91 ± 0.82	88.81 ± 1.47	87.91 ± 0.73	0.76 ± 0.01	0.94 ± 0.00	g = 0.5, c = 2, j = 1
C5 bin	− 0.2	91.18 ± 0.92	86.61 ± 1.68	88.80 ± 1.14	0.78 ± 0.02	0.96 ± 0.00	g = 0.5, c = 1, j = 2
N5C5 bin	0.2	89.20 ± 1.11	91.14 ± 1.61	90.22 ± 1.05	0.80 ± 0.02	0.96 ± 0.00	g = 0.05, c = 1, j = 4
N10 bin	− 0.2	86.39 ± 2.73	89.68 ± 1.79	88.42 ± 1.05	0.76 ± 0.02	0.94 ± 0.01	g = 0.1, c = 2, j = 2
C10 bin	− 0.2	79.87 ± 2.30	86.49 ± 2.43	83.96 ± 1.91	0.66 ± 0.03	0.90 ± 0.01	g = 0.05, c = 3, j = 1
N10C10 bin	− 0.4	86.89 ± 2.70	91.62 ± 2.92	89.83 ± 1.31	0.79 ± 0.02	0.96 ± 0.00	g = 0.1, c = 1, j = 1
AAC + motif	− 0.1	95.51 ± 0.86	95.35 ± 0.85	95.42 ± 0.77	0.91 ± 0.01	0.99 ± 0.00	g = 0.001, c = 6, j = 1
DPC + motif	0	94.15 ± 0.92	96.94 ± 0.49	95.71 ± 0.38	0.91 ± 0.00	0.99 ± 0.00	g = 0.0005, c = 1, j = 2

This table shows average performance (mean ± standard deviation) of models on randomly generated training datasets (bagging)

MCC Matthews correlation coefficient, AAC amino acid composition, DPC dipeptide composition, N5 first 5 residues from N terminus, C5 first 5 residues from C terminus, N5C5 first 5 residues from N and C terminus respectively, bin binary profile, AAC + motif amino acid composition with MERCI motif score, DPC + motif dipeptide composition with MERCI motif score, SVM parameters g gamma parameter of the radial basis function, c trade-off between training error and margin, j regularization parameter (cost-factor, by which training errors on positive examples outweigh errors on negative examples)

Composition-based models

The amino acid composition (AAC) and dipeptide (DPC) composition were used to develop SVM-based models on the training–testing dataset. On evaluating the performance parameters, the AAC model gave an accuracy of 93.30% and the Matthews correlation coefficient (MCC) value of 0.87 as given in Table 1. The DPC model also gave a similar performance in terms of accuracy (94.84%) and MCC (0.90) values. When the AAC and DPC models on the terminal 5 residues of the sequences individually (N or C terminus—N5, C5) or together (N and C termini combined—N5C5) were developed, the N5C5 AAC and N5C5 DPC models performed closest to but not better than the AAC and DPC models with MCC values for N5C5 AAC being 0.87 and N5C5 DPC 0.86 (Table 1).

Binary models

Binary models take the residue position into account by representing each residue type as a binary vector. We considered 5 and 10 residue positions from either end of the peptide sequences (N and C termini) and developed SVM-based models individually and in combination. In case of 5-residue position consideration, the model developed on combined 5 residue positions on both the N and C termini (N5C5 bin in Table 1) performed the best giving an accuracy value 90.22% and MCC value 0.80 while the N10C10 bin model performed the best among the 10-residue position models (accuracy 89.83% and MCC 0.79).

Hybrid model

In the previous motif analysis, motifs exclusive to the A-cell epitopes were identified using the MERCI program. We checked whether the motif information added to composition could help improve the performance of prediction. Indeed, the AAC + motif model that combined the information of presence of motifs exclusive to A-cell epitopes with amino acid composition achieved a better performance than AAC model on the training–testing dataset giving an accuracy value of 95.42% and MCC value of 0.91 (Table 1). Yet, the best model among all the feature combinations was that of motif information combined to the DPC (DPC + motif model) that gave an accuracy of 95.71% and the MCC value 0.91.

Performance of the models on the independent datasets

The independent dataset was generated in the same way as training dataset (by sampling) and consisted of ~ 20% of the total dataset resulting in 61 positive sequences and 77 negative sequences. The average performances of the AAC and DPC models on the independent dataset were comparable to those on the training–testing dataset with AAC model giving an accuracy of 93.91% and MCC of 0.88 while DPC model giving an accuracy of 94.64% and MCC of 0.89 (Table 2). The MCC values of the N5C5 AAC and N5C5 DPC models in the independent dataset evaluation (0.88 and 0.88 respectively) were also close to those found in the training–testing dataset evaluation. Similar to the training–testing dataset results, the N5C5 bin and N10C10 bin models performed the best (MCC 0.82 and 0.81 respectively) among the binary models on the independent dataset. The hybrid model (AAC + motif and DPC + motif) gave accuracies 94.35 and 95.00% respectively, while the MCC values found were 0.89 and 0.90 respectively when evaluated on the independent dataset (Table 2). Figure 5 is a plot of the MCC values of various SVM models on the training–testing dataset along with the MCC values of the models obtained on the independent dataset drawn shown as bars. The MCC values on both the datasets are comparable for each of the models developed indicating that the models are not over optimized on the training–testing dataset.

Table 2

The performance of SVM-based models developed using various features; models were evaluated on independent dataset (external cross-validation)

Features	Threshold	Sensitivity (%)	Specificity (%)	Accuracy (%)	MCC	AUROC	Parameters
AAC*	− 0.1	94.10 ± 2.70	93.77 ± 3.28	93.91 ± 2.00	0.88 ± 0.03	0.98 ± 0.00	g = 0.001, c = 3, j = 3
N5 AAC	0	89.75 ± 3.61	90.88 ± 3.67	90.32 ± 2.27	0.81 ± 0.04	0.95 ± 0.01	g = 0.0005, c = 2, j = 1
C5 AAC	0	91.12 ± 3.53	91.64 ± 2.90	91.40 ± 2.13	0.83 ± 0.04	0.97 ± 0.01	g = 0.001, c = 9, j = 1
N5C5 AAC	− 0.2	94.61 ± 3.35	93.59 ± 3.23	94.07 ± 2.16	0.88 ± 0.04	0.98 ± 0.00	g = 0.0005, c = 1, j = 1
DPC	0	93.77 ± 2.76	95.32 ± 1.40	94.64 ± 1.24	0.89 ± 0.02	0.99 ± 0.00	g = 0.0005, c = 1, j = 2
N5 DPC	− 0.1	81.68 ± 4.25	87.36 ± 2.78	84.62 ± 2.65	0.69 ± 0.05	0.93 ± 0.01	g = 1e−05, c = 9, j = 1
C5 DPC	− 0.1	92.31 ± 3.36	94.71 ± 2.45	93.55 ± 1.40	0.87 ± 0.02	0.98 ± 0.01	g = 0.0005, c = 1, j = 2
N5C5 DPC	− 0.2	94.10 ± 3.06	93.75 ± 2.33	93.90 ± 1.49	0.88 ± 0.03	0.98 ± 0.01	g = 0.0001, c = 1, j = 1
N5 bin	− 0.1	88.46 ± 2.90	89.43 ± 3.32	88.98 ± 2.36	0.78 ± 0.04	0.95 ± 0.01	g = 0.5, c = 2, j = 1
C5 bin	− 0.2	93.70 ± 3.03	87.88 ± 4.25	90.63 ± 2.43	0.82 ± 0.04	0.97 ± 0.01	g = 0.5, c = 1, j = 2
N5C5 bin	0.2	90.95 ± 3.18	91.13 ± 3.44	91.03 ± 2.77	0.82 ± 0.05	0.97 ± 0.01	g = 0.05, c = 1, j = 4
N10 bin	− 0.2	89.38 ± 6.68	90.46 ± 4.67	90.01 ± 3.26	0.79 ± 0.06	0.95 ± 0.03	g = 0.1, c = 2, j = 2
C10 bin	− 0.2	85.02 ± 8.02	85.24 ± 5.15	85.19 ± 4.09	0.69 ± 0.09	0.93 ± 0.03	g = 0.05, c = 3, j = 1
N10C10 bin	− 0.4	88.73 ± 5.95	92.33 ± 5.69	91.04 ± 2.52	0.81 ± 0.05	0.97 ± 0.02	g = 0.1, c = 1, j = 1
AAC + motif	− 0.1	93.11 ± 1.86	95.33 ± 3.13	94.35 ± 1.67	0.89 ± 0.03	0.99 ± 0.00	g = 0.001, c = 6, j = 1
DPC + motif	0	93.28 ± 2.38	96.36 ± 1.70	95.00 ± 1.25	0.90 ± 0.02	0.99 ± 0.00	g = 0.0005, c = 1, j = 2

This table shows average performance (mean ± standard deviation) of models on randomly generated independent datasets (bagging)

MCC Matthews correlation coefficient, AAC amino acid composition, DPC dipeptide composition, N5 first 5 residues from N terminus, C5 first 5 residues from C terminus, N5C5 first 5 residues from N and C terminus respectively, bin binary profile, AAC + motif amino acid composition with MERCI motif score, DPC + motif dipeptide composition with MERCI motif score, SVM parameters g gamma parameter of the radial basis function, , c trade-off between training error and margin, j regularization parameter (cost-factor, by which training errors on positive examples outweigh errors on negative examples)

Fig. 5

Comparison of the support vector machine-based prediction models on the training–testing and the independent datasets. The striped bars correspond to the Matthews correlation coefficient (MCC) values obtained for the models on the training–testing dataset and the solid line joins the MCC values of the models on the independent dataset. For each model, the MCC values for the training–testing dataset and the independent dataset are comparable indicating the reliable prediction capabilities of the models

Alternate negative dataset of random peptides

In the absence of the experimentally verified non-epitopes, the human endogenous circulating peptides were considered to be negative sequences. There are two major issues with this dataset; firstly, it is possible that some of endogenous circulating peptides are A-cell epitopes and secondly, the size of negative dataset is small. In this study, we developed models on alternate dataset. In alternate dataset, non-epitopes were derived randomly from human proteins available in the Swiss-Prot database. Our alternate dataset the non-epitopes (or random peptide) 10 times the number of A-cell epitopes. Bagging procedure of sampling was performed to create ten training and ten independent datasets for internal and external validation respectively. The DPC and DPC + motif models performed better than other models on the main dataset. As shown in Table 3, the performance of models on alternate dataset is similar to the performance of models on main dataset. Overall, the performance values are comparable indicating that the choice of negative dataset in the absence of the experimentally verified negative sequences or the total number of training examples has a minimal effect on the performance of the prediction models developed in this study.

Table 3

The average performance of A-cell epitope prediction models on training and independent dataset

Features	Threshold	Sensitivity (%)	Specificity (%)	Accuracy (%)	MCC	Parameters
Internal validation: performance on training dataset, evaluated using fivefold cross-validation
DPC	− 0.2	87.49 ± 1.41	98.70 ± 0.16	97.68 ± 0.22	0.86 ± 0.01	g: 0.0005, c: 1, j: 4
DPC + Motif	− 0.2	87.81 ± 1.01	99.30 ± 0.10	98.25 ± 0.17	0.89 ± 0.01	g: 0.0005, c: 1, j: 4
External validation: performance on independent dataset
DPC	− 0.2	87.54 ± 4.31	98.87 ± 0.28	97.84 ± 0.41	0.87 ± 0.02	g: 0.0005, c: 1, j: 4
DPC + Motif	− 0.2	77.86 ± 5.84	99.28 ± 0.30	97.33 ± 0.58	0.83 ± 0.04	g: 0.0005, c: 1, j: 4

These training and independent datasets were created from alternate datasets using bagging. In alternate dataset, negative or non-epitopes were derived from human proteins. The performance values have been reported as mean ± standard deviation for each model

MCC Matthews correlation coefficient, DPC dipeptide composition, DPC + motif dipeptide composition with MERCI motif score, SVM parameters g gamma parameter of the radial basis function, c trade-off between training error and margin, j regularization parameter (cost-factor, by which training errors on positive examples outweigh errors on negative examples)

VaxinPAD: the web interface for prediction of A-cell epitopes

For providing the SVM-based A-cell epitope prediction methods developed in the present study to the scientific community, we designed an in silico platform VaxinPAD (available at http://webs.iiitd.edu.in/raghava/vaxinpad/). The platform has a couple of utilities that may help the user in designing peptide-based adjuvants as well enhance or diminish the immunomodulatory potential of the query sequence.

Designing of vaccine adjuvants

The ‘PREDICTION’ module of the VaxinPAD platform allows the user to check whether the query sequence would be immunomodulatory on the basis of SVM score. It allows the virtual screening of the A-cell epitopes among a library of input peptide sequences.

Designing analogs of adjuvant peptides

The ‘ANALOGS’ module of VaxinPAD enables a user to generate all possible single residue position substituents of a query peptide sequence and predict potential immunomodulators among the analogs generated.

Immunomodulatory regions in a protein

The ‘PROTEIN ADJUVANTS’ module of VaxinPAD does a window search across the length of a query protein sequence to identify immunomodulatory patches, the window size being user-defined. LL37, a well-known immunomodulatory peptide is a 37 amino acid peptide derived from human Cathelicidin. This menu may help the researchers in identifying more such peptides that are immunomodulatory.

Peptide sequences dataset

Finally, VaxinPAD includes a menu ‘DATASET’ that includes a list of immunomodulatory peptides collected from literature. Among the sequences in the database only the peptides of length 3–30 were used for development of prediction models in the current work.

Discussion

Previously, peptide-based vaccine adjuvants were largely being developed as ligands of innate immunity receptors like TLR-4 and TLR-2 [45] or as self-assembling nanostructure forming entity [46]. Recently, it has been realized that short immunomodulatory peptides can be developed as potential vaccine adjuvants [4]. Cationic host defense peptides were previously known to have antibacterial activity by direct killing of the pathogen [47]. Of late, these peptides have been found to evoke the innate immunity by a variety of mechanisms [6]. A majority of these mechanisms involve pattern recognition receptors (PRRs) playing important roles especially in the antigen presenting cells (APCs) like dendritic cells, macrophages, etc. Since these peptides activate APCs, we call these peptides as ‘A-cell epitopes’ (antigen presenting cell epitopes). The A-cell epitopes undertaken in the present study were collected from the patent literature that included host defense peptides as immunomodulatory sequences. To the best of the authors’ knowledge, the present study is the first attempt to develop an in silico tool for designing innate immunomodulatory peptides as the first step towards engineering novel peptide-based vaccine adjuvants.

An important finding in this study was that the residues preferred in A-cell epitopes include arginine (R). Arginine enrichment of the peptide sequences is an important aspect of increasing the cellular uptake of cell-penetrating peptides (CPPs) [48]. Cathelicidins are recognized as an important class of host defense peptides that includes many arginine-rich peptides [49]. Further, human cathelicidin-derived peptide LL37 that is rich in basic residues arginine and lysine has been reviewed as a promising immunomodulatory peptide with cell penetrating properties [50]. Hence, sequence analysis of the A-cell epitopes may indicate cell-penetrating ability to be an associated property of the A-cell epitopes.

Another aspect of our sequence analysis is the occurrence of n-mers found sparsely in the naturally occurring proteins. Patel et al. [51] found that peptide pentapeptides occurring rarely in the universal proteome when introduced into the end of the antigenic sequence enhanced its antigenicity and also suggested that on exogenous addition these rare pentapeptides could act as immunomodulators and thus could be developed as adjuvants. In our analysis too, we found tripeptides, tetrapeptides and pentapeptides occurring rarely in Swiss-Prot proteins to be present more in the A-cell epitopes than non-epitopes. This fits well with the intuition that the immune system is more likely to react to rarely encountered sequence motifs than frequent ones.

On evaluating the performance of SVM models based on composition, the dipeptide composition showed no improvement over amino acid composition. The binary models also showed a lower performance than composition-based models. On the other hand, addition of the motif information increased the performance of both the amino acid composition model as well as dipeptide composition to achieve the maximum accuracies of ~ 96%. In our study, we implemented our best models- the dipeptide composition-based and that based on the hybrid of dipeptide composition and motifs. These models were better than other models possibly owing to the fact that the dipeptide composition provides more information in comparison to simple composition. Dipeptide composition provides the information about the amino acid fraction as well as their local order [34]. It is well known fact that there are certain patterns/motifs present in the proteins/peptides which are responsible for its biological activity [52-55]. In past also, many methods have been developed which have shown that adding the motif information increases the accuracy of model [56, 57]. In our analysis, we have identified certain motifs which are exclusively present in A-cell epitopes. We also observed that adding the information of these motifs with dipeptide composition improved the accuracy of the model. Therefore, we believe that this model will help in classifying the A-cell epitope from non-A cell epitope more accurately in comparison to other models.

We also checked whether the random distribution of the main dataset into the dataset used for training the models (internal dataset) and that evaluating them independently (external dataset) renders a bias in the performance of prediction models. Further, the choice of sequences assumed to be the non-epitopes could also affect the performance of the prediction models. A low number of positive sequences in the total dataset could be a third source of influence on the robustness of prediction models. Using various methods, we observed that all of these three factors have a negligible effect on the performance of the best performing SVM-based prediction models for A-cell epitopes.

The peptides designed using the tools developed in the present study might act by various mechanisms and receptors for activating the innate immunity owing to the fact that the training dataset of the prediction models contains peptides acting by diverse cell signaling routes. Hence, the in silico tool presented here could help an investigator to begin with a choice of peptides that may be the starting molecules for the development of vaccine adjuvants. However, there are certain limitations associated with this model, for example, the method does not consider modifications (e.g., post-translational modifications) and other topological aspects during model development. Secondly, whether the predicted peptides actually prove to be useful as adjuvants, would have to be tested experimentally. Another limiting aspect of the present study is the exclusion of very long immunomodulatory sequences and lastly, the size of the dataset used in the study. The model can be further optimized by incorporating more peptide features such as physico–chemical properties, modifications, etc. Also, with larger datasets and receptor-specific ligands made available in future, studies subsequent to the present investigation might help design peptides eliciting a specific desired innate immune response leading to adjuvancy. Nonetheless, VaxinPAD developed in this study for predicting the immunomodulatory peptides sets a stage for the advancement of rational peptide-based vaccine adjuvant designing.

In silico methods for predicting and identifying DNA and RNA-based immunomodulatory molecules have already been developed [13, 14]. The current study is the first attempt to develop models for predicting immunomodulatory peptides for the development of vaccine adjuvants. Though other biomolecules like lipopolysaccharides and glycosaminoglycan’s also cause activation of innate immunity by binding to the PRRs, the literature currently does not hold a sufficient number of molecules for developing prediction models in these categories. Future studies might focus on the development of in silico tools for predicting such immunomodulatory biomolecules for obtaining new vaccine adjuvants. In addition to this, peptides with non-natural chemical modifications might offer better adjuvants too. Correspondingly, computational tools for prediction of modified peptides might also become an area of development.

Conclusion

Host Defense Peptides have been realized as promising immunomodulators likely to become potential vaccine adjuvants [47]. With immunomodulatory properties, novel peptides predicted to be A-cell epitopes using the models developed here are also likely to have potential to provide host protection against pathogens. Many host defense peptides (HDPs) with known immunomodulatory effects are already in clinical trials [47]. Despite the associated toxicity of the A-cell epitopes due to their pleiotropic effects on the immune system, rational design of innate defense regulators (IDRs) that are synthetic analogs of HDPs is in pressing demand for having immunopotentiators with reduced toxicity and increased specificity of immune responses. We have developed SVM-based models for prediction of A-cell epitopes that could be used to formulate vaccine adjuvants. These models have been implemented in the form of webserver VaxinPAD available at http://webs.iiitd.edu.in/raghava/vaxinpad/ and http://crdd.osdd.net/raghava/vaxinpad/ freely to the scientific community.

Additional file 1: Table S1. Dataset (A-cell epitopes and non-epitope) used for training, testing and validation; including source of information. Table S2. The percent amino acid composition and difference in composition for A-cell epitopes, non-epitopes and human proteins including Students t-test with p-value and adjusted p-value. Table S3. The percent dipeptide composition and difference in composition in A-cell epitopes, non-epitopes and human proteins. Table S4. The percent tripeptide compositions of A-cell epitopes, non-epitopes and human proteins. Rows are sorted in decreasing order of values in the 5th column. Table S5. Motifs exclusively found in A-cell epitopes and their frequency of occurrence, sorted in the decreasing order of counts. Table S6. Motifs exclusively found in non-epitopes and their frequency of occurrence, sorted in decreasing order of counts. Table S7. Performance of various classifiers on the training–testing dataset and the independent dataset. Table S8. The performance of SVM-based A-cell epitope prediction models developed using various features; models were evaluating on training and independent datasets using fivefold cross-validation.