Genetic diversity and neutral selection in Plasmodium vivax erythrocyte binding protein correlates with patient antigenicity

Plasmodium vivax is the most widespread and difficult to treat cause of human malaria. The development of vaccines against the blood stages of P. vivax remains a key objective for the control and elimination of vivax malaria. Erythrocyte binding-like (EBL) protein family members such as Duffy binding protein (PvDBP) are of critical importance to erythrocyte invasion and have been the major target for vivax malaria vaccine development. In this study, we focus on another member of EBL protein family, P. vivax erythrocyte binding protein (PvEBP). PvEBP was first identified in Cambodian (C127) field isolates and has subsequently been showed its preferences for binding reticulocytes which is directly inhibited by antibodies. We analysed PvEBP sequence from 316 vivax clinical isolates from eight countries including China (n = 4), Ethiopia (n = 24), Malaysia (n = 53), Myanmar (n = 10), Papua New Guinea (n = 16), Republic of Korea (n = 10), Thailand (n = 174), and Vietnam (n = 25). PvEBP gene exhibited four different phenotypic clusters based on the insertion/deletion (indels) variation. PvEBP-RII (179–479 aa.) showed highest polymorphism similar to other EBL family proteins in various Plasmodium species. Whereas even though PvEBP-RIII-V (480–690 aa.) was the most conserved domain, that showed strong neutral selection pressure for gene purifying with significant population expansion. Antigenicity of both of PvEBP-RII (16.1%) and PvEBP-RIII-V (21.5%) domains were comparatively lower than other P. vivax antigen which expected antigens associated with merozoite invasion. Total IgG recognition level of PvEBP-RII was stronger than PvEBP-RIII-V domain, whereas total IgG inducing level was stronger in PvEBP-RIII-V domain. These results suggest that PvEBP-RII is mainly recognized by natural IgG for innate protection, whereas PvEBP-RIII-V stimulates IgG production activity by B-cell for acquired immunity. Overall, the low antigenicity of both regions in patients with vivax malaria likely reflects genetic polymorphism for strong positive selection in PvEBP-RII and purifying selection in PvEBP-RIII-V domain. These observations pose challenging questions to the selection of EBP and point out the importance of immune pressure and polymorphism required for inclusion of PvEBP as a vaccine candidate.


Introduction
In 2017 Plasmodium vivax caused between 7.5-14.3 million cases of malaria, mainly in the South-East Asia region (56%) [1,2]. Although P. vivax is generally not lethal to their host, P. vivax causes high morbidity among the five human invasive Plasmodium species (P. falciparum, P. vivax, P. knowlesi, P. malariae and P. ovale) due to recurrent parasitaemia from reactivation of the dormant hypnozoites [1,[3][4][5]. Despite its importance as the most widespread Plasmodium species, absence of reliable in vitro long-term culture system has impeded research into optimal interventions against the parasite [5,6]. Current malaria vaccine development strategies focus on identifying a specific, immunogenic antigen which will stimulate protective humoral immune response and to produce sufficient amount of the specific, functional antibody to provide sterile immunity against malaria infection. Finding such functional antibodies from individuals living in endemic settings may answer for natural infections against malaria. Anti-malarial humoral immune responses provide various function including phagocytosis and/or direct killing by complement mediation [7] which results in a reduction in merozoite invasion, growth and rosetting formation [7][8][9].
Erythrocyte invasion of Plasmodium species occurs by sequential multiple molecule interactions, with each step mediated by antigens belonging to different protein families [10]. The erythrocyte binding-like (EBL) protein family which contained EBL (or RII) domain was identified in various Plasmodium species with conserved function as host cell binding via glycoprotein receptor for invasion before tight junction formation [11,12]. These proteins are highly expressed in mature schizont stage and localized in the microneme [12]. The best-known member of EBL family is P. vivax Duffy binding protein (PvDBP) which binds to Duffy antigen receptor for chemokine (DARC/Fy glycoprotein) and mediate invasion into erythrocytes by the majority of P. vivax isolates [13]. This specific ligand-receptor interaction and parasite invasion process can be a target for invasion blocking antibodies against EBL domain (PvDBP-RII) [14,15]. Given its unique importance, the EBL domain has been the main candidate for vaccine development. However, recent Phase 1a PvDBP-RII vaccine clinical trial showed low efficacy [16] and strain-specific immune response is limited by a high level of genetic polymorphism in the EBL domain [17]. To overcome this important issue, identifying conserved epitopes and evaluating their recognition by neutralising antibodies is essential to achieve sterile immunity [18,19]. Further investigations of novel antigens for vaccine candidates are required to understand downstream immune response of target antigen, evaluation of current immunity to conserved epitopes, and determination of regions undergoing neutral selection.
Here, we have focused on another member of EBL family, i.e. P. vivax erythrocyte binding protein (PvEBP), which was first identified from a Cambodian field isolate (C127) [20]. PvEBP preferentially binds CD71 + and CD234 + reticulocytes but not CD234reticulocyte through its EBL domain, PvEBP-RII [21]. Antibodies against PvEBP-RII but not PvDBP is inhibited binding of reticulocytes to PvEBP expressing COS cells indicating that despite sequence similarity the binding regions of these two EBL proteins are different [21]. An increasing number of studies have reported P. vivax invasion in Duffy-negative populations in Africa regions highlighting potential of PvDBP independent invasion pathways [22][23][24]. Given its similarity, it is hypothesised that, PvEBP may have a role in facilitating this newly discovered DARCindependent invasion pathway [21]. In addition to the analysis of PvEBP for its functional contribution to parasitic invasion, there is also a need to evaluate immune responses and natural polymorphisms of EBP in field isolates to further determine the likelihood of EBP as a possible vaccine candidate.
We analysed clinical isolates from eight countries to quantify genetic diversity of PvEBP and total serum concentrations of IgG against two domains of PvEBP (RII and RIII-V). Our analysis evaluates the correlation of vivax patient antigenicity with genetic polymorphism and natural selection to explore the prevalence of PvEBP specific IgG and how this might impact vaccine development.

Nucleotide diversity and neutral selection
Pvebp nucleotide diversity (π) is defined as the average number of nucleotide differences per site between two sequences within the sequences. The number of polymorphic sites, number of haplotypes and haplotype diversity (Hd) were calculated by DnaSP software [28]. Test of neutral selection was evaluated using multiple calculation method including Tajima's D, Fu and Li's D � , and Fu and Li's F � with excluding the sequence gap [29,30]. Under neutrality, Tajima's D is expected to be 0. Significant positive values of Tajima's D indicate population bottlenecks and balancing selection, whereas negative values suggest population expansion or negative selection [29]. Significant positive value of Fu and Li's D � and Fu and Li's F � represent population contraction due to selection. On the other hand, negative values represent population expansion and excess of singletons [30]. Natural selection was determined by calculating the rates of synonymous substitutions per synonymous site (d S ) and non-synonymous substitutions per nonsynonymous site (d N ) at the intra-species level. The calculation was performed by computed Nei and Gojobori's method and robustness was estimated by the bootstrap method with 1000 pseudo replicates as implemented in the MEGA 5 software. A d N /d S ratio less than 1 indicates a purifying selection and d N /d S ratio greater than 1, indicates a positive selection. The test of pvebp natural selection at the inter-species level was performed using the robust McDonald and Kreitman (MK) test with P. cynomolgi DBP2 (PcyM_0102400) gene as an out-group using DnaSP software. PvEBP haplotype was evaluated by DnaSP software and graphical presentation for distance in relationship was generated by median-joining method in Network 5.0 software.

Antigenicity evaluation
Protein microarray was performed to evaluate total IgG reactivity. 3 aminopropyl-coated slides were prepared as described previously [31]. The slides were printed to each spot with recombinant protein (RII, 200 ng/μl and RIII-V, 12.5 ng/μl) as a saturated concentration and incubated for 2 hours at 37˚C. The recombinant protein coated slide was blocked with blocking buffer (5% BSA in PBS-T) for 1 hour at 37˚C. Vivax patient and healthy individual sera were diluted in PBS-T to 1:25 and probed on the chip for 1 hour at 37˚C. The arrays were visualized with 10 ng/μl of Alexa Fluor 546 goat anti-human IgG (Invitrogen, Carlsbad, CA, USA) in PBS-T for 1 h at 37˚C and scanned with ScanArray Gx laser confocal scanner (PerkinElmer, Norwalk, CT, USA). The positive cut-off values calculated by negative control mean fluorescence intensity (MFI) plus two standard deviations.

Statistical analysis
The data were analyzed using GraphPad Prism (GraphPad Software, San Diego, CA, USA), SigmaPlot (Systat Software Inc., San Jose, CA, USA), and Microsoft Excel 2013 (Microsoft, Redmond, WA, USA). For the protein array, Student's t-test was used to compare the experimentally measured values of each group. The correlation of clinical information with antigenicity was calculated by Pearson correlation test. Differences of p < 0.05 were considered significant. The total IgG reactivity index was calculated by each MFI divided by average negative MFI of RII and RIII-V, respectively, for normalization and reactivity comparison between RII and RIII-V.
Pvebp sequences identified 67 haplotypes with four distinct clusters. The cluster 1 contained 61 haplotypes with geographically sharing for eight countries (Fig 2). The result indicates that the pvebp cluster 1 polymorphism occur widely without specific geographical pattern. The far distance from core haplotypes population were revealed specific geographical pattern from Ethiopia for cluster 2 (H51), cluster 3 (H46, H50, and H51), and ROK isolates for cluster 4 (H01 and H02) (Fig 2).

Antigenicity screening
PvEBP-RII (30.1 kDa) and PvEBP-RIII-V (25.8 kDa) recombinant proteins were used for antigenicity evaluation (Fig 4A). Total IgG reactivity in vivax patient showed significantly higher in RII (MFI ± S.D., 16,048 ± 8,547) than RIII-V (7,386 ± 3,211). Additionally, RII IgG reactivity in healthy individual also showed significantly higher in RII (11,175 ± 5,010) than RIII-V (5,462 ± 1,762) (Table 3). However, PvEBP specific IgG inducing level after P. vivax infection  in RIII-V (21.5%) was higher than RII (16.1%) ( Table 3). This result indicated that RII antigen mainly recognized by natural IgG for innate protection, whereas specific IgG inducing level by B-cell was lower than RIII-V domain. On the other hand, RIII-V had lower level of natural IgG reactivity caused by neutral selection for immune evasion, whereas PvEBP specific IgG production level was higher than RII after infection ( Fig 4B and Table 3). Antigenicity correlation analysis between RII and RIII-V in individual vivax patient showed significant correlation from ROK (r = 0.710), Thailand (r = 0.534), and Myanmar (r = 0.305) (Fig 4C). Combination of all countries correlation between RII and RIII-V showed low level of significant positive correlation (r = 0.417) (Fig 4). The correlation analysis of RII and RIII-V with age and parasitaemia showed positively correlation only in RIII-V with patient age (Fig 5A-5D).

Discussion
This study describes the genetic diversity of the extracellular domain (ecto) of P. vivax erythrocyte binding protein (pvebp) genes from eight countries and correlated these polymorphisms with vivax patient antigenicity from different endemic areas. The genetic phenotypes of pvebpecto were classified into four clusters, based on the insertion/deletion (indels) variation. Generally, the gene indels variation of pathogen directly affects host antibody recognition. Cluster 1 is a major genetic phenotype of pvebp which widely shares its geographical haplotype. However, three minor clusters (cluster 2, 3, and 4) which provide geographically distinct patterns. Interestingly, the large insertion sequence in PvEBP-RIII-V domain of cluster 3 and 4 was conserved in all of the tested field isolates regardless of geographical location. This indels variation is interesting to note that PvEBP cluster 4 phenotype is closely related with P. cynomolgi Mulligan strain DBP2 (PcyM DBP2) which has a Gly-Lys (GK) amino acid deletion in EBP-RII domain following the large insert in EBP-RIII-V domain. Thus, PcyM DBP2 gene is clearly an orthologue of PvEBP. Recently, P. cynomolgi, has been used as a model to study the invasion biology of P. vivax [34,35]. The EBL family in P. vivax has two genes named PvDBP and PvEBP, whereas P. cynomolgi has three genes named PcyDBP1, PcyDBP2, and PcyEBP (PcyM DBP2) [36]. The additional Gly-and Ser-rich sequence insertion in both P. vivax and P. cynomolgi EBP genes may have a function as a linker that possibly affect flexibility for inducing an action radius of ecto domain and hamper host antibodies recognition [37,38]. However, the function of this insertion on PvEBP and PcyEBP (PcyM DBP2) needs to be determined.
The importance of EBL family proteins RIII-V domain as a vaccine candidate was suggested in the previous study of antibodies against P. falciparum erythrocyte binding antigen 175 (PfEBA-175) RIII-V domain which was associated with protection against symptomatic case   of malaria [39]. PfEBAs (PfEBA-175, PfEBA-140, and PfEBA-181) RIII-V domain was dominantly eliciting IgG1 and IgG3 by patient age-dependent manner [39]. The function of IgG subtype in immune strategy of malaria infection is well documented, in the opsonisation capacity of IgG1 and IgG3 [8]. Similarly, an earlier study demonstrated that an IgG1 response against PvEBP domain (161-641 aa.) was dominant in an age-dependent manner in vivax patients [40]. Our study confirmed that only the IgG response against PvEBP-RIII-V (480-676 aa.) domain was positively correlated with patient age. Thus even though PvEBP-RIII-V domain has large sequence insertion and shows low frequency in the tested isolates, the RIII-V domain needs to determine its function and considered for vaccine design in the future.
The EBL (RII) domain in Plasmodium species is the key domain for host cell binding and invasion by interacting with host cell glycoprotein receptor [11]. Therefore, this functional domain of PvEBP-RII is the most logical target for a parasite-erythrocyte invasion blocking vaccine [41]. Antibodies against PvEBP-RII were shown to prevent the binding of a recombinant protein encoding this domain to reticulocytes [21]. Additionally, antibodies against PvEBP-RII were useful for serological marker of recent malaria exposure, with a half-life estimated at 734 days [42]. Anti-PvEBP-RII antibodies were predominant IgG1 and IgG3 which is similar to other PfEBAs-RII to mediate opsonic phagocytosis [39,40]. However, PvEBP-RII showed highest genetic polymorphism within pvebp-ecto domains. This high polymorphic pattern is consistent with other EBL members in P. vivax and P. falciparum [11,41] and its function as binding domains has been conserved in the other members of the EBL family [43]. As high polymorphisms, PvDBP-RII domain presented difficulty in developing an efficacious vaccine [44,45]. It would be important to consider the genetic polymorphisms of PvEBP-RII when developing a vaccines targeting this domain [17]. Study of PvEBP copy number variation (CNV) revealed higher multiple copy ratios in Madagascar (56%) than in Cambodian (19%) isolates, possibly correlating with the Duffy dependence phenotype [46]. Geographically, PvEBP gene CNV pattern was similar to PvDBP which detected higher CNV in Africa region such as Ethiopia (79%) and Madagascar (46%) than South-East Asia such as Thailand (30%), Cambodia (28%), Indonesia (6%), and Malaysia (4%) [26,46]. In this study, Ethiopia isolates covered three different clusters (clusters 1, 2, and 3) which indicate high level of genetic and phenotypic diversity in the African region where the population is largely Duffy-negative [47]. The Duffy negativity may have affected both pvebp and pvdbp gene-phenotype to change according to host cell preference or/and immune evasion mechanism [23]. The neutrality test of pvebp-ecto showed that rare alleles were present at high frequencies. Especially, PvEBP-RIII-V which showed evidence of population expansion and a negative/purifying selection effect. This was in contrast to the PvEBP-RII domain which showed positive selection pressure. Previous reports of PvEBP-RII neutral selection revealed positive/diversifying selection similar to other EBL domains in human and non-human primate malaria [46,48,49]. All of these results support that adaptation to the host environment lead both PvEBP-RII and PvEBP-RIII-V domains to be highly variable.
The antigenicity screening result also reflects the influence of a high level of a PvEBP gene variant. Although both PvEBP-RII (16.1%) and PvEBP-RIII-V (21.5%) domains significantly induced antibody responses in vivax infected patients, these responses were comparatively lower than other P. vivax antigens [50][51][52]. Previous studies on vivax malaria patients, (using the same antigenicity screening method used in this study) conveyed the antibody responses to be PvDBP-RII (56.9%), PvRBP1a-34 (33.7%), PvRBP1b-32 (39.4%), and PvGAMA-ecto (72.0%) which were localized at apical organelle and important for host cell interaction [51,52]. Relatively low antigenicity of PvEBP could be a result of high polymorphism with positive selection and population expansion in the PvEBP-RII domain, and purifying selection in PvEBP-RIII-V domain. Interestingly, total IgG production and recognition of PvEBP-RII and PvEBP-RIII-V showed clearly distinct properties. The total IgG response against PvEBP-RII was higher than PvEBP-RIII-V, however, acquired antibody response was less than PvEBP-RIII-V. In contrast, PvEBP-RIII-V elicited higher IgG responses in vivax patients, yet, basal IgG recognition level was poor.
PvEBP is a prominent novel vaccine candidate for blood-stage P. vivax malaria with its potential to target both Duffy-dependent and independent P. vivax parasites. However, low antigenicity due to high genetic polymorphism in the PvEBP-RII, and antigen phenotype and selection pressure in the PvEBP-RIII-V will need to be considered for future vaccine development.