Rubella virus is a positive sense single strand RNA virus that belongs to Rubivirus genus of the Togaviridae family, it contains 9.7 KB genome which surrounded by capsid protein (nucleocapsid protein). The nucleocapsids surrounded by envelop that contains E1 and E2 glycoproteins, these proteins arranged in parallel rows on the outer surface, giving the virion different shapes [1-3]. Rubella virus is classified into two clades, clade 1 and clade 2 according to E1 protein, and they are divided further into 11 genotype (1B, 1C, 1D, 1E, 1F, 1G, 1H, 1I, 1J, 2A, 2B) and 1 provisional genotype [4-6]. After outbreak that occurred in 1964 in USA and Europe which affected about 1% of pregnancies, rubella infection was studied more [7]. Worldwide there are more than 100000 babies affected by congenital rubella syndrome each year, and there were more than 800000 rubella cases reports in 1999 and 2001 with majority of cases reported by Russian federation [8]. Rubella is highly infectious disease transmitted by respiratory droplet, usually causes mild disease in children with symptoms of low grade fever, maculopapular rash, arthralgia and myalgia [9-11]. However infection of women during pregnancy especially in first 12 weeks may cause death of fetus or Congenital Rubella Syndrome which has different clinical manifestations such as sensory neural deafness, congenital cataract and glaucoma, microcephaly and heart disease that occur in 80%-90% of cases [9,10].The immune response in human body developed after vaccination and infection, with formation of antibodies against structural protein [12-15], which directed mainly toward E1, E2 and capsid protein. The antibodies against E1 glycoprotein  persist for long time, in contrast, the antibodies against E2 and capsid protein tend to persist for short time, and this is may be due to the difference in accessibility by  immune system [2,5,16-19].

The current live attenuated vaccine is Rubella vaccine RA 27/3 (human diploid fibroblast) strain which used in combined with mumps and measles vaccine (MMR) and effectively reduce disease incidence [20-23]. However, the vaccine has some side effects such as joint symptoms that may occur in up to 60% of post pubertal female but does not cause chronic joint disease, also lymphadenopathy, rash and allergy [24-26]. In fact, it is contraindicated in immunosuppression, malignant disease, and during pregnancy because of its risk to the fetus, but there was no confirmation of congenital rubella syndrome in babies of pregnant women who took the vaccine [10,27-29]. Because of significant effects of rubella on people health, a safe vaccine with maximum efficacy, good coverage and least side effects is needed, especially for immunosuppressed patients and pregnant women due to high disease burden. Thus the aim of this study was to determine antigenic peptides from E1, E2 and Capsid proteins that can be used for multiple peptide based vaccine design against rubella virus using immunoinformatics analysis. This study is unique in rubella virus vaccine design because it deals with three structural proteins that can induce protective immune response, in order to cover most possibilities of epitopes that can be used as vaccine to achieves maximum results, Unlike other type of vaccine, using In-Silico studies to identify epitopes will insure less risk of virulence and side effects that can be seen with lives attenuated vaccine [30], as only specific epitope will be used, it also allow production of vaccine with chemical safety, low cost and good population coverage as only conserved epitopes are tested.

A total of 183 rubella virus E1 glycoprotein sequences, four E2 glycoprotein sequences and only two capsid protein sequences were obtained from NCBI (www.ncbi.nlm.nih.gov/protein/?term=) in FASTA format in March 2017 (Table 1).

Table 1: Countries and accession numbers of retrieved sequences from NCBI;
*Ref sequence, remaining data as extra file.

For all three proteins, reference sequences were identified from NCBI reference sequence as shown in Table 2.

Table 2: Reference sequence of three proteins; AA: Amino Acid.

All sequences of each protein were subjected to multiple sequence alignments using CLUSTALW tool of BIOEDIT sequence alignment editor (version 7.2.5.0) in order to identify conserved regions between sequences (Figure 1). Then epitopes prediction and analysis of each protein was done using different tools of immune epitope data base IEDP software (http://www.iedb.org) [31].

Figure 1: Multiple sequence alignment of E1 protein sequences from data base using BIOEDIT software. Uncolored regions represent mutated residues in the sequence.

Epitopes Prediction

Prediction of B cell epitopes: The first step in the identification of antigenic epitopes in the pathogens is the identification of linear peptide parts. A combination between hidden Markov model (HMM) and (Parker and Levitt) method was done to predict epitopes accurately (Figure 2) [32]. Using IEDB software particularly Bepipred linear epitope prediction tool, B cell epitopes from conserved regions were identified from the E1, E2 and capsid protein with specific default threshold value for each protein, as shown in Table 3 [33].

Table 3: Default threshold for all three proteins used in Linear B cell epitopes prediction, surface accessibility and epitopes antigenicity prediction.

Figure 2: Prediction of B cell epitopes by various scales for E1 protein; Remaining figures as extra file.

Surface accessibility prediction: The indices of surface probability method were developed to increase confidence in provisional alignment for comparing other sequences predicted by the chou and fasman method and garnier et al. method, And this surface probability method assumes the absence of significant internal deletions or insertions [34]. Emini surface accessibility prediction tool of IEDP was used to predict surface accessibility using default threshold value for each protein as shown in Table 3 [34].

Epitopes antigenicity sites: Identification of antigenic sites for each protein with default threshold value was conducted using Kolaskar and Tongaonker antigenicity tool of IEDB as demonstrate in Table 3 [35].

MHC epitope prediction: IEDP server (http://www.iedb.org) was used through specific tools to determine MHC1 and MHC II binding epitopes. This server uses specific scoring IC50 (inhibitory concentration 50) to predict epitopes that bind to different MHC class I and MHC class II alleles.

MHC class I epitope prediction: The MHC system is an example of receptor that can interact with linear ligands of variable lengths, and the method of prediction for MHC class I affinity was been tested on large set of quantitative peptide MHC class I measurement affinity on the IEDB [36], by using artificial neutral network (ANN) method and length of nine amino acid, all conserved epitopes bound with score equal or less than 300 IC50 for all three structural proteins were chosen for further analysis [37].

MHC class II binding prediction: One of the major goals of the immunological research is the identification of MHC class II restricted peptide epitopes, and for that reason many computational tools have been developed but also their performance has lack of large scale systematic evaluation, and we use comprehensive dataset consisted from thousands of previously unpublished MHC peptide binding affinities and peptide MHC crystal structures, which all tested for CD4+ T cell responses to evaluate the publicly available MHC class II binding prediction tools performances [38]. Therefore MHCII binding tool from IEDP was used by applying NN align as prediction method. Then using IC50 prediction value equal or less than 1000, all conserved epitopes were chosen for more analysis [38].

Population coverage: The estimation of the population coverage are based on the MHC binding with or without T cell restriction data, there for Nemours based tool was developed to predict population coverage of T cell epitope-based diagnostic and vaccines based on MHC binding with or without T cell restriction data [39]. All alleles that interact with epitopes from E1, E2 and capsid protein were subjected population coverage tool of IEDB (http://tools.iedb.org/tools/population/iedb_input) to calculate the whole world population coverage of MHC class I, MHC II and combined MHC I and II alleles for each protein [39].

Homology modeling: The 3D structure for the three different proteins was obtained using Raptor X structure prediction server [40] and Chimera 1.8 [41] was used to demonstrate the structure of proposed B cell and T cell epitopes that can be utilized for vaccine development as shown in Figure 3.

Figure 3: Proposed epitope for vaccine development; a) Proposed peptide of E1 protein that bind to B cell (magenta), MHC I (red), MHC II (white); b) Proposed peptide of E2 protein that bind to B cell (blue), MHC I (magenta), MHC II (magenta& white); c) Proposed peptide of capsid protein that bind to B cell (white), MHC I (magenta), MHC II (blue).

Prediction of B cell epitopes

All three proteins (E1, E2 and capsid protein) was tested using Bepipred linear epitope production, Emini surface accessibility and Kolaskar and Tongaonkar antigenicity tool of IEDB. In Bepiprd linear epitope prediction method the average binding score for E1, E2 and capsid protein were 0.32, 0.26 and 0.97 respectively, with minimum and maximum values as shown in Table 4. All values equal or greater than the default threshold were considered as potential B-cell binders.

Table 4: Average, minimum and maximum score values of the three tests.

Regarding Emini surface accessibility prediction, the average binding score of E1, E2 and capsid was 1.00, and minimum and maximum values as displayed in Table 4. All values equal or greater than the default threshold were predicted to have good surface accessibility.

All proteins were subjected to Kolaskar and Tongaonkar antigenicity prediction tool of IEDB to predict peptides with probability of being antigenic, the average threshold value for E1 was 1.05, E2 was 1.04 and 0.90 for capsid protein, with minimum and maximum values shown in Table 4. Epitopes with values equal or greater than the average score are considered as antigenic peptides.

For E1 protein seven predicted epitopes had succeeded both antigenicity and surface accessibility test. Peptide PVCQRHSP from 233 to 240 was found to have high score.

The result of E2 glycoprotein showed that seven peptides had passed the antigenicity prediction and surface accessibility prediction test, AQYPP from 163 to 167 was found to have greatest score among all other peptides in both tests.

Analysis of Capsid protein resulted in nine peptides that passed both tests, and PPHT from 258 to 261 was found to have high score in Emini surface accessibility and Kolaskar and Tongaonkar antigenicity prediction test. The result is summarized in Table 5. Various tools illustrated in Figure 2.

Table 5: Result of predicted B cell epitopes by different tools; Remaining data as extra file.

Prediction of T-cell epitopes and MHC class I interaction analysis

Epitopes sequences of E1, E2 and capsid subjected to MHC class I binding prediction tool of IEDP. T-cell epitope predicted to interact with different MHC class I alleles using ANN (artificial natural network) as prediction method and length of nine amino acids. Using same score for all protein of 300 IC50, 12 peptide of E1 were found to interact with MHC class I. The peptide YFNPGGSSY (101-109) had height affinity to interact with the largest number of alleles (10 alleles, HLA-A*29:02, HLA-A*30:02, HLA-B*15:01, HLA-B*15:02, HLA-B*35:01, HLA-B*35:01, HLA-C*03:03, HLA-C*07:02, HLA-C*12:03, HLA-C*14:02) when compared to other peptides.

For E2 protein 52 peptides were predicted to interact with MHC class I alleles, the peptide RRRGAAAAL from 260-268 was found to have high affinity to interact with largest number of alleles (HLA-B*07:02, HLA-B*27:05, HLA-B*39:01, HLA-C*07:02, HLA-C*14:02) followed by the peptide 242 FVLLVPWVL 250 which predicted to interact with three alleles (HLA-A*02:01, HLA-A*02:06, HLA-B*39:01).

Analyses of capsid protein resulted in 54 peptides interact with different MHC class I alleles. Among these the peptide FTNLGTPPL from 174-182 was interacted with highest number of alleles (8 alleles, HLA-A*02:01, HLA-A*02:06, HLA-A*68:02, HLA-B*35:01, HLA-B*39:01, HLA-C*03:03, HLA-C*05:01, HLA-C*14:02). All data are summarized in Table 6.

Table 6: Predicted peptides that interacted with MHC class I alleles from E1, E2 and capsid protein; Remaining data as extra file.

Prediction of T Helper cell epitopes and interaction with MHC class II

Using MHC class II binding prediction method of IEDB and based on NN align and IC50 of 1000 for E1, E2 and capsid protein. From all predicted epitopes of E1 protein that bind to MHC class II, the core sequence FTGVVYGTH bind to 11 alleles (HLA-DPA1*01:03, HLA-DQA1*05:01, HLA-DRB1*04:01, HLA-DRB1*04:04, HLA-DRB1* 04:05, HLA-DRB1*09:01, HLA-DRB1*11:01, HLA-DRB5*01:01, HLA-DPA1*01:03, HLA-DPB1*02:01, HLA-DQB1*03:01) which is largest number of alleles in contrast to other predicted core sequences.

From Analysis of E2 protein interaction with MHC class II, the core FVLLVPWVL was found to bind to 16 alleles (HLA-DRB1*01:01, HLA-DRB1*04:05, HLA-DRB1*07:01, HLA-DRB1*15:01, HLA-DQA1* 05:01, HLA-DRB1*01:01, HLA-DRB5*01:01, HLA-DPA1*01, HLA-DPA1*01:031, HLA-DPA1*01:03, HLA-DPA1*03:01, HLA-DQB1* 02:01, HLA-DPB1*02:01, HLA-DPB1*02:0, HLA-DPB1*04:02, HLA-DPB1*04:01), while prediction of capsid protein binding to MHC class II revealed that the core FYRVDLHFT bound to the highest number of 15 different alleles (HLA-DRB1*01:01, HLA-DRB1*04:01, HLA-DRB1*04:05, HLA-DRB1*07:01, HLA-DRB1*11:01, HLA-DPA1* 01, HLA-DPA1*01:03, HLA-DPA1*02:01, HLA-DRB1*01:01, HLA-DPA1*03:01, HLA-DPB1*04:01, HLA-DPB1*04:02, HLA-DPB1* 05:01, HLA-DPB1*01:01, HLA-DPB1*02:01). The result is summarized in Table 7.

Table 7: Result of predicted peptides from structural proteins that interact with MHC class II alleles; Remaining data as extra file.

Analysis of population coverage

All epitopes from the three structural proteins which selected as epitopes with high affinity to interact with MHC class I and MHC class II, subjected to IEDB population coverage set against the whole world population.

For MHC class I, Epitopes with highest population coverage were YFNPGGSWW (54.08%), FVLLVPWVL (42.23%) and FTNLGTPPL (58.29%) for E1, E2 and capsid protein respectively. Data displayed in Table 8.

Table 8: Population coverage of proposed peptide interacted with MHC class I from the three proteins; Remaining data as extra file.

In MHC class II, the epitopes that interact with class II alleles and showed highest population coverage were found to be FTGVVYGTH (93.05%), FVLLVPWVL (99.95%) and FYRVDLHFT (99.93%) for E1, E2 and capsid respectively, as shown in Table 9. The epitope FVLLVPWVL form E2 protein is interacted with both MHC I and MHC II with high affinity and highest population coverage for each class. All these proposed peptides are illustrated in Figure 3 at structural level.

Table 9: MHC class II Population coverage of proposed peptides; Remaining data as extra file.

The combined coverage of MHC I and MHC class II

Using IEDB and population coverage prediction tools, the top three peptides that interact with most frequent alleles for both MHC class I and class II for each protein were analyzed and the result showed epitope set of 98.54% and 99.99% for E1 and E2 respectively, while combined population coverage of epitopes from capsid protein displayed an epitope set of 100%. The data demonstrated in Table 10.

Table 10: Combined population coverage of MHC class I and class II for proposed peptides from each protein.

Prediction of population coverage was also calculated for combined MHC class I and class II using proposed epitopes from ALL three protein (E1, E2 and capsid protein) which resulted in epitope set of 100%. The data was displayed in Table 11.

Table 11: Combined population coverage of MHC class I and class II for proposed peptides from all three proteins E1, E2 and capsid protein.

The aim of this study is to identify potential epitopes that can induce cellular and humoral immune reaction and act as candidate for rubella vaccine development. We used Immuno informatics tools to identify epitopes for multiple peptide vaccine, these tools was recently used to design vaccines for many viruses such as Ebola [42], Zika virus [43], merkel cell polyomavirus [44] and human papillomavirus [45]. Rubella virus contains two projections that are mainly protein in nature, E1 and E2, which anchored to the outer layer of the membrane and many forms of nucleocapsid envelope protein. The E1 protein is the immunodominant antigen and also plays role in endocytosis and viral neutralization [46-48]. Babies with congenital rubella syndrome showed proliferative response against E2, and significant proportion showed proliferative response against capsid protein, while antibodies response was mainly directed against E1 followed by E2 and capsid protein [49,50]. Immunoglobulin has been used in attempts to prevent rubella in pregnant women exposed to the virus. However, it does not appear to be highly effective [46]. The first step in the infection is the membrane fusion of E1 an E2 with the enhancement of capsid protein [51]. Fusion of the viral envelope occurs when conformational change of E1 and E2 proteins induced by exposure to pH of 6 or less [52]. Our proposed peptides form three proteins are 100% conserved regions in proteins as well as succeeded all subjected tests for B cell and also bound to MHC alleles with high affinity.

Rubella E1 protein was subjected to IEDB B cell epitope prediction tests. It was found that the most satisfactory peptide is 8 amino acid PVCQRHSP B cell epitope from 233 to 240 with antigenicity score of 1.116 and 1.167 Score for Emini surface accessibility. Mitchell et al. measured specific IgG antibodies level before and after reimmunization with MMR vaccine and found similar result as E1(234-252) peptide contain antibody neutralizing domain [12], while E1(208-239) bind effectively to monoclonal antibodies in study done by Wolinsky et al. using mice and rabbit immunized by this peptide [53]. While 12 Epitopes from E1 Protein interacted with MHC class I HLA alleles. The proposed T cell peptide 101 YFNPGGSSY 109 is well conserved among Rubella E1 protein. It was identified using IEDB MHC I prediction tool. Yang et al. found that the peptide between the 81 and 109 in E1 protein involved in the membrane fusion activity of the virus [54] reflected the importance of this protein for the virus. Chong et al. reported that T cell epitopes are mainly located in E1 and capsid protein but less in E2, but they recognized different peptides from E1, E1(207- 226), E1(324-343) and E1(358-377) that can stimulate cellular immune response [55], our peptide predicted to bind to 10 different HLA alleles (HLA-A*29:02, HLA-A*30:02, HLA-B*15:01, HLA-B*15:02, HLA-B*35:01, HLA-B*35:01, HLA-C*03:03, HLA-C*07:02, HLA-C*12:03, HLA-C*14:02) with high affinity by MHC prediction tool. Ovsyannikova et al. analyzed HLA class I type in group of children who received second dose of MMR and they reported other alleles, their result showed that HLA-B 3503 and HLA-CW 1502 were associated with T cell response against the vaccine [56]. All epitopes were tested using population coverage tool of IEDB which measure the percentage of people in whole world who have potential to develop immune response to vaccine contains this epitope. Our proposed peptide that binds to MHC 1 alleles had coverage of 54.08% which is the highest figure among all predicted peptides. On another handout of seven Predicted peptides our proposed core FTGVVYGTH is part of different peptide sequence in E1 protein from 417 to 434, interact with 11 MHC class II HLA alleles with high affinity, Mitchell et al. found increase in simulation indices after MMR vaccine for different E1 peptides, E1(213-239), (234-252), (254-285), (272-285), (301-314) and (462-481) [12]. Ovsyannikova et al. tested the association between HLA class II and cellular and humoral immune response after rubella vaccine in 346 children and found that DPB1*0301, DQB1*0501, DRB1*0101, and DRB1*1104 associated with cellular immune response [57], which are different from what we predicted using this epitope. Our core sequence has potential population coverage of 93.05% which considered as promising coverage for vaccine that will contain this epitope. Rubella peptides induce protective neutralizing antibodies [46], that considered to be protective in contrast to measles vaccine which does not prevent infection or disease [58].

Among all predicted E2 epitopes from only conserved region as in Table 5 which subjected to Bepipred linear epitope prediction, Emini surface accessibility and Kolaskar and Tongaonker antigenicity, the epitope AQYPP from 163 to 167 Had highest score in both test of 2.363 in surface accessibility and 1.074 in antigenicity test, and chosen as proposed peptide from E2 protein that can activate B cell to produce antibodies against the virus. In contrast, Mitchell et al. reported that E2 (1-16) and E2 (10-36) had domain that can neutralize antibodies [12]. T cell immune response is essential for longer lasting response [59]. For E2 protein, among all 52 conserved predicted epitopes which interact with high affinity to MHC class I as summarized in Table 6, the epitope FVLLVPWVL from 242 to250 was found to had high affinity to interact with three MHC class I alleles, which is dissimilar to E2(1-16), (10-36), (35-58), (50-72), (134-150), (140-156), (168-179) and (248-260) that stated to have cellular immune response as reported by Mitchell et al. [12]. Our proposed peptide has highest world population coverage of 42.23%. The same epitope had also the highest affinity to interact to different 16 MHC class II alleles (HLA-DRB1*01:01, HLA-DRB1*04:05, HLA-DRB1*07:01, HLA-DRB1*15:01, HLA-DQA1*05:01, HLA-DRB1* 01:01, HLA-DRB5*01:01, HLA-DPA1*01, HLA-DPA1*01:031, HLA-DPA1*01:03, HLA-DPA1*03:01, HLA-DQB1*02:01, HLA-DPB1* 02:01, HLA-DPB1*02:0, HLA-DPB1*04:02, HLA-DPB1*04:01), out of 59 predicted core sequences. Ovsyannikova et al. found association between HLA class II alleles and T cell response specifically DRB1 (DRB1*0101, DRB1*0701, and DRB1*1104), DQB1(DQB1*0202 and DQB1*0501) and DPB1 (DPB1*0301,DPB1*0401, DPB1*1001 and DPB1*1101) response to MMR vaccine [57] which have some alleles in common with our study. This core resulted in population coverage of 99.59% as predicted by IEDB and 99.97% for combined MHC I and MHC II worldwide coverage. The greatest results for both classes make it promising epitope for vaccine design as it able to induce T cell immune response.

The capsid protein is part of host interaction and important for virus assembly [51]. It undergoes a structural change permitting release of viral genomic RNA into the cytoplasm [52]. Nine epitopes from only conserved region was found to interact with B cell, epitope258PPHT261 represented the proposed epitope with high Emini surface accessibility and antigenicity score of 1.285 and 1.036 respectively. Few studies conducted about capsid protein as protein for vaccine design, but Lovett et al. reported that capsid peptide C (1-29) contain T cell epitope that help in antibodies production [60]. In comparing to all 54 epitope from conserved regions which interacted with MHC I alleles, epitope 174 FTNCGTPPL 182 bind to 8 alleles (HLA-A*02:01, HLA-A*02:06, HLA-A*68:02, HLA-B*35:01, HLA-B*39:01, HLA-C*03:03, HLA-C*05:01, HLA-C*14:02). Chong et al. recognized T cell epitopes in capsid (119-152), (205-233) and (255- 280) [55] which similar to some predicted peptides but different from our proposed peptide. With potential population coverage of 58.29%, this peptide considered as candidate for vaccine production. The core FYRVDLHFT has high score among 51 conserved epitopes that interacted with MHC II alleles by binding to 15 Alleles. In contrast to Capsid peptide (11-29) that was promiscuously recognize HLA class II restricted CD4 T cell as reported by Lovett et al. [60] our peptide resulted in very high population coverage of 99.95%. To get maximum benefits of vaccine we predicted the population coverage of combined T cell epitopes in Tables 6 and 7 for each protein and the results were promising, with epitope set of 98.54% for E1, 99.99% for E2, and 100% for capsid protein. Furthermore vaccine that contains the top three epitopes from all structural protein predicted to have 100% coverage in the whole world reflecting the favorable effect of this vaccine.

The efficacy and safety of predicted epitopes by this computational analysis are needed to be evaluated by animal model studies, to confirm whether they can induce protective immune response or not. There are only six sequences of E2 and capsid proteins available in the database; more sequences are needed to increase the significance of the result. The following proposed peptides are recommended for multiple peptides vaccine design against rubella virus; (PVCQRHSP from 233 to 240, YFNPGGSYY from 101 to 109 and the core FTGVVYGTH) from E1. (AQYPP from 163 to 167 and FVLLVPWVL from 242 to 250) from E2, (PPHT from 258 to 261, FTNLGTPPL from 174 to 182, and the core FYRVDLHFT) from capsid protein. This vaccine will insure good population coverage and fewer side effects that can be seen with life attenuated vaccine.

The efficacy and safety of predicted epitopes by this computational analysis are needed to evaluate animal model studies, to confirm whether they can induce protective immune response or not. There are only six sequences of E2 and capsid proteins available in the database; more sequences are needed to increase the significance of the result. The following proposed peptides are recommended for multiple peptides vaccine design against rubella virus; (PVCQRHSP from 233 to 240, YFNPGGSYY from 101 to 109 and the core FTGVVYGTH) from E1. (AQYPP from 163 to 167 and FVLLVPWVL from 242 to 250) from E2, (PPHT from 258 to 261, FTNLGTPPL from 174 to 182, and the core FYRVDLHFT) from capsid protein. This vaccine will insure good population. Coverage and fewer side effects that can be seen with life attenuated vaccine.