Plasmodium vivax HAP2 genetic diversity and population structure from worldwide clinical samples. A potential Transmission-Blocking malaria vaccine candidate

1 Introduction

In 2020, an estimated 241 million malaria cases were reported worldwide, with 6.2 million deaths (World Health Organisation 2022). Plasmodium vivax is a significant public health problem, disproportionately impacting the poorest and most vulnerable populations. The overall proportion of cases has reduced from 8 % in 2000 to 2 % in 2020, according to the World Malaria Report (2022). Investments in malaria control and elimination programs have significantly reduced global malaria cases in recent decades. However, in areas where both P. falciparum and P. vivax coexist, P. vivax is increasingly becoming the primary cause of malaria (Battle et al., 2019). Even though P. vivax is widely distributed, there are presently no effective vaccinations available (Antonelli et al., 2020). Of all documented cases, 29–75 % of P. vivax malaria occurs occurs in the regions of the Americas, Southeast Asia, and the Eastern Mediterranean area (Poespoprodjo et al., 2009). The parasite has a complex, multi-stage life cycle, which makes developing effective vaccines for Plasmodium species difficult. For example, a wide and stage-specific protein repertoire is expressed during the asexual and sexual life cycle stages of P. vivax parasites, intermediate hosts, namely mosquitoes, and definitive hosts, which are humans (Florens et al., 2002; Le Roch et al., 2004). It is therefore necessary to develop innovative and integrated malaria prevention and control measures. A highly efficacious malaria vaccine is crucial to control and eventually eradicate malaria. Transmission-blocking vaccines (TBVs) that decreases the transmission of the parasite from humans to mosquitoes represent a significant step towards the objective of eliminating malaria. TBVs stimulate the production of antibodies that specifically target malaria antigens found in sexual stages, including gametocytes, gametes, zygotes, and ookinetes. These antibodies play a vital role in reducing the density of oocysts in the mosquito's midgut. In particular, transmission-blocking vaccinations (TBVs) are regarded as a critical element in efforts to control and eradicate malaria worldwide (PARASITE, 2002; Saul, 2007; Wu et al., 2015; Acquah et al., 2019). A TBV may be an effective technique for the eradication of P. vivax infections because of the early and continuous generation of gametocytes during infection. This is a crucial factor for the transmission of the parasites via mosquitoes. TBVs target antigens in the malaria parasite's sexual and early sporogonic phases, as well as midgut proteins in the mosquito vectors (Wu et al., 2015; Acquah et al., 2019). Proteins which are involved in fertilization are identified as targets for the development of TBVs and several of these have been well characterised in P. falciparum which shows transmission-reducing activity and these include P48/45, P47, P25, P28, P230, AgAPN1 and gamete membrane fusion protein HAP2/Generative Cell specific 1 (GCS1) (Duffy and Kaslow, 1997; van Dijk et al., 2001; Dinglasan et al., 2007; Wu et al., 2008; Blagborough and Sinden, 2009; Marin-Mogollon et al., 2018). Among these, HAP2/c, P230 and P48/45 are pre-fertilization antigens, while P25 and P28 are post-fertilization antigens and AgAPN1 is a mosquito midgut antigen.

HAP2/GCS1 is categorized as a class II viral fusion protein and contains a cysteine-rich. domain. Initially, this protein was discovered in the plant species Arabidopsis thaliana and Lilium longiflorum (Mori et al., 2006; von Besser et al., 2006; Angrisano et al., 2017a). The HAP2/GCS1 domain is found in a wide variety of organisms, including Plasmodium species and in the green algae Chlamydomonas. In Plasmodium, HAP2/GCS1 protein is localized on the surface of male gametes (Mori et al., 2006; Liu et al., 2008). A conserved 'cd loop' sequence of the HAP2/GCS1 protein, which corresponds to amino acids (a.a.) 178–207 in P. falciparum and (a.a.) 174–205 in PbHAP2, is essential for membrane binding during fertilisation (Angrisano et al., 2017b; Fédry et al., 2017). Mouse antibodies against recombinant PfHAP2 expressed in a wheat germ cell-free system showed strong transmission-reducing activity (Miura et al., 2013). A recent study reported that anti-PvHAP2 antibodies substantially reduced the number of oocysts in mosquitoes suggesting that PvHap2 could be a potential TBV against P. vivax infections and aid in malaria elimination (Qiu et al., 2020). Altogether, these studies suggest that HAP2 is a potential transmission-blocking vaccine candidate however, very limited studies have been done to investigate genetic diversity studies from worldwide samples. A recent study from the Greater Mekong Subregion (GMS) conducted on a fragment from 13nt to 2574nt of the gene indicated limited genetic diversity among field samples within the region (Li et al., 2020a). Similarly, a study focusing only on the short GCSI domain i.e. 397nt – 2097nt of the HAP2 gene from Iran indicated limited genetic diversity too (Mehrizi et al., 2016). Therefore, in this study, the genetic diversity, natural selection, haplotype network analysis and population structure of PvHAP2/GCS1 full-length gene using worldwide samples from 16 countries are determined and the implications for its use as a potential vaccine candidate against P. vivax are discussed. Potential epitopes within different domains of the protein are identified too.

2 Materials and Methods

3 Results

3.1

3.1 PvHAP2/GCS1 sequence identity, polymorphisms, and phylogenetic relationship between p. Vivax HAP2/GCS1 and its ortholog species

Alignment and comparison of deduced amino acid sequences of P. vivax HAP2/GCS1 (Sal I strain) with its ortholog members in P. knowlesi, P. falciparum, P. reichenowi, P. ovale and P. malariae was found to be 85.7 %, 51.1 %, 51.1 %, 62.1 % and 64.9 % respectively. A schematic diagram of the full-length structure of the HAP2/GCS1 gene of P. vivax and 30 amino acid polymorphisms with their precise locations within each domain as observed within the 118 worldwide sequences are shown in (Fig. 1). The study revealed 18 PvHAP2 amino acid haplotypes from the 118 worldwide amino acid sequences, haplotype 3 (n = 35, 29.6 %) was shared among 11 countries (Columbia, Mexico, Peru, Thailand, India, China, Gabon, Madagascar, Myanmar, PNG, and Cambodia) (Figure S1). Other shared haplotype with lower number of samples included haplotype 1 (n = 25, 21 %), haplotype 4 (n = 24, 20 %) and haplotype 2 (n = 20, 16.9 %) (Figure S1). The phylogenetic analysis using deduced amino acid sequences indicated that PvHAP2/GCS1 is closer to its ortholog in the zoonotic malaria parasite P. knowlesi compared to its human infecting orthologs in P. falciparum, P. malariae and other primate malaria species i.e. P. coatneyi (Fig. 2).

Close

Fig. 1

Domain organisation of PvHAP2 and amino acid polymorphism observed within 118 worldwide PvHAP2 sequences obtained from 16 countries. SP; Signal peptide, DI, DII, DIII, the domain I-III; L, linker (L); TM, transmembrane segment; C-rich region, the cysteine-rich region; PFAM10699, the HAP2/GSC1 domain.

Export to PPT

Close

Fig. 2

Phylogenetic relationship of PvHAP2 with its orthologs in Plasmodium knowlesi, Plasmodium malariae, Plasmodium coatneyi, Plasmodium falciparum and Plasmodium reichenowi. The silhouette figures of primate and humans indicates the host for each of the parasites.

Export to PPT

3.2

3.2 PvHAP2 diversity and natural selection

DNA sequence analysis of 118 PvHap2 sequences revealed that there were 95 single nucleotide polymorphisms (SNPs) within 118 PvHAP2 sequences from 16 different nations, of which 65 were synonymous substitutions and 30 were non-synonymous substitutions (Figure S1). The hydrophobic “fusion loop” (cd loop) segment (167–198 a.a.) did not contain any non-synonymous substitution indicating conservation of membrane fusion function and gamete fertilisation. The overall nucleotide diversity (within all countries) was found to be π = 0.0030 which was slightly higher compared to a previous study done with samples only from the Greater Mekong Sub region (Li et al., 2020a) (Table 1). There were 68 parsimony informative sites within PvHAP2, 27 were singleton variable sites, and there were 25 haplotypes with a haplotype diversity of 0.842. (Table 1). The graphical representation of the nucleotide diversity and Tajima’ D values of PvHAP2/GCS1 is shown in (Fig. 3A) and (Fig. 3B) respectively. Thirty non-synonymous substitutions identified within the amino acid sequence are shown in (Figure S1).

Close

Fig. 3

(A) Graphical representation of nucleotide diversity (π) within PvHAP2 gene and (B) Tajima’s D value. The window length and step size of π and D graph is 100 & 25 respectively as implemented in DnaSP software v5.0.

Export to PPT

Natural selection analysis of PvHAP2/GCS1 gene from 118 sequences indicated the gene is under negative or purifying selection (dN-dS = -5.10, P > 0.1) probably due to functional constraints (Table 1). Overall, the Tajima’s D values were negative and significant (D = − 1.79, P < 0.05) indicating a strong purifying selection. The Fu & Li D* and F* values were also negative but not significant which indicates purifying selection and population expansion (Table 1).

Table 1 Genetic diversity and neutrality indices of PvHap2 sequences from worldwide sequences.

Domain	No. of Samples	SNPs	Syn	NonSyn	No. of haplotypes	Diversity ± SD		Codon based z-test dN─Ds	Fu & Li’s D*	Fu & Li’s F*	Taj D
						Haplotype	Nucleotide
Full length	118	95	65	30	25	0.842 ± 0.018	0.0030 ± 0.00078	−5.1 (p > 0.1)	−1.78 (p > 0.10)	−1.24 (p > 0.10)	−1.79 (p < 0.05)

3.3

3.3 Genealogical haplotype network analysis of PvHAP2/GCS1

The nucleic acid based haplotype network analysis of 25 haplotypes of the PvHAP2/GCS1 revealed that there were 4 major haplotype clusters (H_2, H_3, H_4 and H_9), of which, H_2, H_3 and H_4 formed a closely related group compared to Hap_9 indicating the presence of two major haplotype groups (Fig. 4). Among the major haplotypes, H_2 was shared within 11 countries and it had the highest number of samples (n = 34) (Fig. 4). The haplotype 9 (H_9) however, formed a separate cluster with maximum samples from Myanmar, China and Gabon. Minor shared haplotypes (H_7 and H_6) were observed within major haplotype clusters (H_2, H_3, H_4 and H_9).

Close

Fig. 4

Haplotype network of PvHAP2 from 118 sequences from 16 countries. The highlighted clusters represents the two different subpopulation. Various colors represent different locations. Circles represent haplotypes with sizes proportional to haplotype frequency.

Export to PPT

3.4

3.4 Population structure analysis of PvHAP2/GCS1 gene from 16 countries

A Bayesian algorithm-based population structure analysis of 118 PvHAP2/GCS1 sequences from 16 countries using STRUCTURE software revealed that there were two distinct sub-populations of the parasite within the worldwide samples (K = 2, ΔK = 503.64) (Fig. 5A and B).

Close

Fig. 5

Most likely number of Plasmodium vivax parasite subpopulation clusters (K) within 16 countries based on Hap2 sequences (K = 2, ΔK = 503.64). Bayesian model–based STRUCTURE version 2.3.4 software was used.

Export to PPT

3.5

3.5 B cell epitopes identified in PvHAP2

Based on amino acid haplotype analysis, the full-length PvHAP2 sequence of amino acid haplotype 3 was shared among 11 countries (FigureS1) for determining the B cell epitopes. Bcpred and IEDB servers identified 19 and 16 epitopes in P. vivax HAP2 amino acid sequence, respectively (Table 2) (Fig. 6). Epitope outputs from both servers were compared, and the results revealed that there are a total of 11 PvHAP2 epitopes, with their lengths varying from 5 to 26 peptides identified by both servers (Table 2). Among the common peptides identified, a lysine-rich repeat peptide was identified which had a high score from both the servers (782 KQTRRKKKKKKS 796).

TABLE 2 The possible epitope predicted by using IEDB server and Bcrepd server using P. vivax HAP2 amino acid sequence is shown in table. The common peptides within each epitopes observed are bolded.

Bcpred server				IEDB server
No.	Start AA	End AA	Peptide	Length	Start AA	End AA	Peptide	Length
1	1	5	MKPTR	5	21	27	RDKREDG	7
2	19	29	FGRDKREDGGR	11	117	131	HQKFKKYSEREIKEY	15
3	38	48	HSFAKKKVCTS	11	170	180	LFRDNKDIKRS	11
4	97	105	YLKYMKDIP	9	230	236	QEYSYDD	7
5	116	132	NHQKFKKYSEREIKEYT	17	249	256	TKEKKYEL	8
6	172	184	RDNKDIKRSKLKC	13	387	392	QRYYDA	6
7	246	258	DLVTKEKKYELND	13	401	406	GASKYK	6
8	402	410	ASKYKIKNV	9	413	425	SEPQTRIYKPSTL	13
9	415	422	PQTRIYKP	8	431	441	DKINNNHKKID	11
10	426	444	PDHLKDKINNNHKKIDAND	19	622	637	SVTKYDSEPEASQNDD	16
11	505	512	SNQSESKD	8	660	665	YYTKRI	6
12	537	543	TDKETGK	7	729	736	RGKRSKCK	8
13	575	581	NKEIKTM	7	774	780	RVKKSAK	7
14	587	598	RDSNGKECSKED	12	782	796	EALKQTRRKKKKKKS	15
15	624	631	TKYDSEPE	8	808	814	HSRNDYS	7
16	660	672	YYTKRIKNVIKKF	13	838	854	SQSDRSSNKGTRAPSRA	17
17	723	738	KERLHIRGKRSKCKKG	16
18	774	799	RVKKSAKSEALKQTRRKKKKKKSGRM	26
19	840	846	SDRSSNK	7

Close

Fig. 6

The graph represents the possible epitope region in PvHAP2 protein using epitope prediction software IEDB. The peaks (in yellow) in graph shows the possible epitope region.

Export to PPT

4 Discussion and Conclusion

Although the World Health Organisation (WHO) considers that vaccines are one of the most effective strategies for the elimination and eradication of malaria, previous studies have shown that developing malarial vaccines is an extraordinarily difficult task. Two of the most significant obstacles in the process of designing a vaccine are allele-specific immune responses (Fluck et al., 2004) and the antigenic diversity of vaccine candidate antigens (Girard et al., 2007; Arnott et al., 2012). To design and develop a malaria vaccine that is effective on a global scale, it is vital to investigate the genetic diversity of vaccine candidate antigens in both low-endemic and high-endemic regions. The surface antigens of gametes and ookinetes are the primary focus of research and development for TBVs. In this study, the genetic diversity, phylogeography, population structure, and epitope characterization of Plasmodium vivax HAP2 full-length gene sequences were determined from worldwide clinical samples (n = 118, 16 countries) obtained from public databases. Previous population genetic studies have only focused on the GCSI domain (Mehrizi et al., 2016; Li et al., 2020a) of the HAP2 protein, which is essential for binding to host cells. Because of this, in the current study, the PvHAP2 full-length gene study from 16 countries was included. Amino acid sequence analysis indicated that the highest similarity was with the monkey malaria parasite P, knowlesi (∼85 %), this might positive implication towards development of cross-species vaccine candidates as reported in previous studies. Though 30 non-synonymous substitutions across the full-length protein were identified interestingly, 15 of these non-synonymous substitutions were localised towards the C-terminal region (674a.a. − 862a.a.) indicating immune selection pressure could be high in this domain. The localisation or accumulation of these non-synonymous substitutions in the C-terminal region could also indicate that this region might be exposed to host immune pressure. However, structural as well as functional studies along with field-based immunological studies need to be conducted to determine the immunogenicity of this domain. Previous studies have focused only on characterizing the HAP2 GCS1 domain of this protein (Qiu et al., 2020). Additionally, genetic diversity studies with a higher number of samples from each location would be essential too.

Analysis of PvHAP2 genetic diversity indices indicated lower levels of diversity (∼0.003) which was slightly higher than previously reported studies (Mehrizi et al., 2016; Li et al., 2020a). This may be due to the fact that this research concentrated primarily on the shorter segments rather than the full-length gene. Other TBVs such as Pvs28 and Pvs25 had shown similar levels (∼0.001––0.004) of genetic diversity in endemic areas of Myanmar and China (Feng et al., 2011; Le et al., 2019) indicating conserved protein function. It is interesting to note that despite low levels of nucleotide diversity, the haplotype diversity was found to be on the higher side for the PvHAP2 gene, indicating population expansion supported by negative values of Li and Fu’ D* and F* (−1.75 & −1.24). Natural selection analysis of worldwide PvHAP2 indicated purifying selection pressure acting on the whole gene which was not statistically significant but was similar to the previous studies and other TBVs (Feng et al., 2011; Mehrizi et al., 2016; Li et al., 2020b).

The haplotype network analysis indicated the presence of four major shared haplotypes within the 16 countries, which were grouped into two major groups. The robust Bayesian model-based population structure analysis also suggested that the PvHAP2 gene largely consisted of two sub-populations (K = 2, K = 503.6). This indicates that any vaccine design needs to take into consideration these two sub-populations to achieve high efficacy and worldwide coverage.

B cell epitope prediction using two different algorithms showed that 11 epitopes were predicted and among these a high scored C-terminal lysine-rich repeat region with 15 peptides EALKQTRRKKKKKKS could be a potential region where high immunogenicity could be found. The lysine-rich 15 amino acid peptide identified in this study could be potential TBV candidate for P. vivax, however, higher number of samples and functional studies would necessary to substantiate the presence of this epitope. This needs to be confirmed with functional as well as field-based immunological studies in the future. Moreover, it should be noted that only two countries, Gabon and Cameroon, both from Central Africa, out of the total 15 countries, have amino acid changes in the C-terminal region (as highlighted in the attached supplementary figure S1). This means that the majority of the field isolates function remained intact, except for parasite from these two countries. We believe that more field-based studies are needed from these countries to confirm whether natural selection is active. In addition, functional studies including sequences from these countries would be essential to understand the function, as we have noted that there are at least 8 amino acid changes in the C-terminal region, which could be important. This is because, the C-terminal region may prove to be crucial for vaccine and therapeutic strategy development, as previous studies on the functional analysis of the cytoplasmic domain of AMA1 have highlighted its significance in invasion and the development of novel therapeutic strategies (Treeck et al., 2009). The lysine-rich 15 amino acid peptide identified in this study could be potential TBV candidate for P. vivax, however, higher number of samples and functional studies would necessary to substantiate the presence of this epitope. Further confirmation through functional and field-based immunological studies on the C-terminal region will be necessary in the future.

5 Conclusion

This work is the first of its kind to disclose the genetic diversity as well as the population structure of the full-length PvHAP2 gene from 17 different nations. Although a low degree of genetic diversity and purifying selection were discovered, the C-terminal region indicated a higher number of non-synonymous substitutions. This finding suggests that this protein is vital to the survival and transmission of the parasite. The haplotype network and population structure analyses both pointed to the existence of at least two distinct subpopulations. This study uncovered a prevalent and shared PvHAP2 haplotype found all around the world, which is notable. In addition, the current investigation was able to identify probable B-cell epitopes over the entirety of the gene, as well as one high-scoring epitope towards the C-terminal region, which consisted of a lysine-rich 15 amino acid peptide. The findings of the current study will be very important in the rational design of a vaccine that can inhibit PvHAP2 transmission and has global application. However, functional studies focusing on immunogenicity and vaccine efficacy studies would be necessary.