Recent increase of genetic diversity in Plasmodium vivax population in the Republic of Korea

Background The reemergence of Plasmodium vivax in South Korea since 1993 represents a serious public health concern. Despite the importance in understanding genetic diversity for control strategies, however, studies remain inconclusive with the general premise that due to low rate of malaria transmission, there is generally low genetic diversity with very few strains involved. In this study, the genetic diversity and population structure of P. vivax in South Korea were explored by analysing microsatellite polymorphism. Methods Sequences for 13 microsatellite loci distributed across the twelve chromosomes of P. vivax were obtained from 58 South Korean isolates collected during two sampling periods, namely 1997-2000 and 2007. The sequences were used for the analysis of expected heterozygosity and multilocus genotype diversity. Population structure was evaluated using STRUCTURE version 2.3.2. Linkage disequilibrium was also analysed to investigate the extent of outbreeding in the P. vivax population. Results Mean expected heterozygosity significantly increased from 0.382 in 1997-2000 to 0.545 in 2007 (P < 0.05). The number of multilocus genotypes was 7 and 27; and genotype diversity was statistically significant (P < 0.01) at 0.661 and 0.995 in 1997-2000 and 2007, respectively. Analysis by STRUCTURE showed a more complex population structure in 2007 than in 1997-2000. Linkage disequilibrium between 13 microsatellites, although significant in both time points, was notably lower in 2007. Conclusions The present microsatellite analysis clearly showed recent increase of genetic diversity and recent relaxation of the strong population structure observed in 1997-2000. These results suggest that multiple genotypes not present previously recently migrated into South Korea, accompanied by substantial outbreeding between different genotypes.


Background
Of the five human malaria parasites, Plasmodium vivax is the most prevalent in Asia, Melanesia, the Middle East, South and Central America, accounting for 70-80 million cases annually [1], with 2.6 billion people at risk of infection [2]. Despite modest gains in Plasmodium falciparum control, the global burden of P. vivax remains underestimated [3] and emergence of drug resistant P. vivax makes the control of vivax malaria more difficult than before [4,5]. In the Republic of Korea (South Korea), vivax malaria was eradicated once in the late 1970s [6]. However, indigenous malaria has reemerged since the first infection of a soldier who had never been abroad was confirmed in 1993 in the Demilitarized Zone (DMZ) between South and North Korea [7]. Thereafter, vivax malaria infection has spread to civilians and the number of cases from 2000 has fluctuated between 864 and 4,142 [8]. Thus, vivax malaria has become endemic again and constitutes a serious public concern in South Korea. Meanwhile, in North Korea, vivax malaria is more prevalent than in South Korea and the number of cases has fluctuated between 7,436 and 296,540 from 2000 [6,9].
Genetic diversity and population structure of P. vivax have a significant impact on malaria transmission, spread of drug-resistance and the acquisition of protective immunity against malaria. Studies of the parasite population diversity have practical significance for the strategic development and deployment of control measures [10]. Limited studies on antigen-encoding genes such as merozoite surface protein-1 (MSP-1), circumsporozoite surface protein (CSP) and merozoite surface protein-3α (MSP-3α) in South Korea have previously shown that the diversity of those genes was low in 1990s, but increased after 2000 [9,11]. An increase in antigen diversity potentially affects the acquisition of protective immunity against vivax malaria. It is therefore important to determine whether the observed recent increase of the diversity of antigenencoding genes is caused by human immune pressure or not. Since these antigens are exposed to the human immune system, diversity in the antigen-encoding loci can be attributed to not only to natural selection by immune pressure but also to the parasite's population history [12]. Hence, it remains uncertain whether the low diversity observed in the antigen-encoding loci resulted from selective pressure imposed by the host immune system or due to a demographic change in the local parasite population. In contrast to polymorphic antigen genes, microsatellite markers are selectively neutral and highly polymorphic, and, thus, are suitable for the assessment of genetic diversity resulting from demography at the genome-wide level. Recently, microsatellite markers have been widely used to analyse genetic diversity of P. vivax populations [13][14][15][16][17][18]. In the present study, multilocus microsatellite analysis was conducted to elucidate genetic diversity and population structure of P. vivax in South Korea. Results showed drastic changes in genetic diversity and population structure between 2000 and 2007.

Parasite samples
Blood samples were collected from indigenous P. vivaxinfected patients between 1996 and 2007 near the DMZ in the Republic of Korea [9]. In this study, two groups classified according to the sampling year were analysed: one group of 29 samples collected between 1997 and 2000, and the other from 29 samples collected in 2007. Informed consent was obtained from all patients before blood sampling. The study was approved by Korea's National Institute of Health, and Center for Disease Control and Prevention.

Preparation of parasite DNA
Plasmodium vivax DNA were extracted from 200 μL blood using proteinase K digestion technique and followed either by purification with phenol-chloroform or QIAamp ® DNA Mini Kit (Qiagen) following the manufacturer's instructions.

PCR and sequencing analysis for genotyping
Thirteen validated polymorphic microsatellite markers, 1.501, 3.27, 3.502, MS3, MS15, MS5, MS9, MS16, MS20, MS6, MS8, MS10 and 14.297 [15,19], distributed across the twelve chromosomes of P. vivax (Additional file 1) were used for analysis. Microsatellites were sequenced to achieve higher-resolution analysis rather than the widely employed fragment-length measurements. Besides, P. vivax microsatellites used in this study contain relatively long and somewhat degenerating repeating units, which cannot be satisfactorily assessed by fragment-length measurements. For example, in MS20 two alleles observed with the same length, (GAA) 5 (CAA) 3 and (GAA) 6 (CAA) 2 , would be indistinguishable from each other by fragment-length measurement. Primers for each locus were newly designed in order to perform polymerase chain reaction (PCR) and subsequent sequencing (Additional file 1). PCR reactions were performed in 20 μL reaction mixture containing 1× LA PCR™ Buffer, 2.5 mM MgCl 2 , 0.4 mM dNTP, 0.2 μM of each primer solution, 0.05 U LA Taq ® (Takara Bio Inc.) and 1 μL DNA template. The reaction protocol was as follows: 1 cycle at 93°C for 1 min; 40 cycles at 93°C for 20 s, 62°C for 4 min; and 1 cycle at 72°C for 10 min. The nucleotide sequence of each locus was determined directly from the PCR products using BigDye Terminator v3.1 Cycle Sequencing Kit (Applied Biosystems) on an Applied Biosystems 3130 Genetic Analyzer (Applied Biosystems) following the manufacturer's instructions. An allele for each locus corresponds to the number of repeat-motifs. Genotype is defined as the unique combination of alleles across the 13 microsatellite loci examined for this study.

Data analysis
Expected heterozygosity (h E ) of each locus was calculated by using the formula h E = {n/(n-1)}(1-∑pi 2 ), where n is the number of isolates analysed and pi is the frequency of the i th allele in the population. Variance (V) of h E was calculated by using the formula V = 2/{n(n-1)}[a+b], where a = 2(n-2){∑pi 3 -(∑pi 2 ) 2 } and b = ∑pi 2 -(∑pi 2 ) 2 [20]. The standard error (SE) of h E was defined as the square root of V. Mean expected heterozygosity (H E ) of all loci and the standard error were also calculated. The statistical significance of the mean differences between expected heterozygosity values was determined by t-test as described by Nei, 1987 [21]. Genotype diversity (gd) was calculated similar to expected heterozygosity.
The standardized index of association (I S A ) was calculated to test the presence of overall multilocus linkage disequilibrium [22]. This test compares the variance (V D ) of the number of alleles shared between all pairs of haplotypes observed in the population (D) with the variance expected under random association of alleles (V E ) as follows: where r is the number of loci analysed. V E is derived from 100,000 simulated data sets in which alleles were randomly reshuffled among haplotypes. Linkage disequilibrium (LD) is significant if V D is > 95% of the values derived from the reshuffled data sets. Data were analysed using LIAN version 3.5 software [23]. LD for all pairs of loci was also tested using the Genepop version 4.1 [24] under the following Markov chain parameters: dememorization number = 20,000, number of batches = 500 and number of iterations per batch = 10,000.
Population structure was evaluated using the STRUC-TURE version 2.3.2 software [25]. This Bayesian modelbased clustering method assigns samples to K populations according to allele frequencies of each locus. The program was run 20 times at each of 10 different K values (1-10) with a period of 50,000 burn-in steps followed by 100,000 iterations. The admixture model with independent allele frequencies was used in all analyses. The number of populations was inferred by plotting the log probability of the data [Ln P(D)] for each K value. Also, in order to determine the most appropriate number of K, ΔK was calculated as described by Evanno et al [26]. Illustrations were generated using the DIS-TRUCT program [27].

Genetic diversity
Microsatellite sequences from 29 isolates in 1997-2000 (Group 1) and another 29 in 2007 (Group 2) were obtained in this study (n = 58 isolates in total). There were two cases in the 1997-2000 group where three sequences could not be obtained. All microsatellite loci were polymorphic, exhibiting three to six alleles per locus (Table 1 (Table 1). All the loci, except 1.501 and MS9, showed higher h E in 2007 than in 1997-2000. These exceptions underscore the importance of using several markers to characterize genetic diversity of a population. Mean expected heterozygosity (H E ) was 0.382 ± 0.061 and 0.545 ± 0.063 in 1997-2000 and 2007, respectively ( Table 1). The difference between the two sampling periods was statistically significant (P = 0.023), indicating a recent increase in the genetic diversity of the P. vivax population in South Korea.
A total of 33 multilocus genotypes, which were defined as the unique association of alleles in the 13 microsatellite loci, were detected: seven in 1997-2000 and 27 in 2007 with one genotype (no. 25) found in both groups ( Table 2). Two of the seven genotypes (no. 25 and 30) comprised 81.5% (22/27) in 1997-2000. Remarkably, for 2007, most genotypes were found only once, except for no. 7 and 8. These suggest a drastic change in genetic diversity between the two groups. Genotype diversity (gd) was 0.661 ± 0.005 for 1997-2000, and 0.995 ± 0.000 for 2007 (P = 0.00004). Hence, the diversity of multilocus genotypes appears to have significantly increased between 2000 and 2007.

Linkage disequilibrium
The overall multilocus linkage disequilibrium assessed by LIAN showed significant LD in 1997-2000 and 2007.

Discussion
This study is the first to report microsatellite diversity of the P. vivax population in South Korea. The present results obtained with 13 microsatellite markers across 12 chromosomes revealed a higher genetic diversity in 2007 than in 1997-2000. Consistent with previous limited studies showing a recent increase in antigen diversity, the study further support changes in P. vivax isolates from South Korea. Analysis of South Korean isolates before 2000 using antigen-encoding genes such as Duffy-binding protein [28], PvCSP [9], PvMSP-3α [11], and AMA-1 [29] infers only two alleles whereas genetic diversity of these antigens remarkably increased after 2001 [9,11]. The present study shows only two major multilocus microsatellite genotypes that comprise the two major P.  [14]).
Further investigations using more samples should be performed to confirm a low genetic diversity in P. vivax population in South Korea. Two mechanisms can be proposed to account for the observed increase in the genetic diversity of P. vivax population in South Korea. One is the accumulation of microsatellite mutations in the population, and the other is migration of novel parasite populations from other geographic areas into South Korea. Since the mutation rate(s) of P. vivax microsatellite is unknown, accumulation of mutation cannot be completely excluded. However, the present study does reveal the appearance of very many novel alleles at almost all the microsatellite loci examined. Previous studies have also reported rapid increase in genetic diversity in antigenencoding genes [9,11]. It is very unlikely that in a short period (between 2000 and 2007, probably six generations of the parasite life cycle in South Korea), many novel alleles accumulate simultaneously at several microsatellite and antigen loci. The idea of migration of multiple different parasite genotypes from other geographic areas into South Korea is, therefore, favoured. Haplotype network analysis based on the parasite mitochondrial genome infers a genealogical origin in southern China for P. vivax populations collected in 1999 and one in 2002 [30]. Further investigations on genetic diversity of P. vivax in neighbouring geographic areas can elucidate and verify the population mechanism.
The rate of inbreeding also contributes to the extent of genetic diversity in P. vivax populations. In this study, the extent of inbreeding inferred from linkage disequilibrium was reduced in 2007, compared with that in 1997-2000; although overall linkage disequilibrium was significant in both sampling groups. These results indicate very strong linkage disequilibrium, leading to a high inbreeding rate in the parasite population in 1997-2000, whereas relatively weak linkage disequilibrium, leading to elevation of outbreeding levels in 2007. Analysis by STRUCTURE ( Figure 2C and 2D) supports the increased outbreeding in 2007. It is thus concluded that the recent diversification of South Korean P. vivax population resulted partly from an increase in the level of outbreeding among different parasite genotypes.
Moving forward, molecular tools and techniques can elucidate parasite genetic diversity and population structure that can help malaria control efforts and understand mechanisms of pathogenicity and drug resistance. It is generally believed in a parasite population with a high genetic diversity, protective immunity against malaria is slow to develop. There is accumulating evidence, which suggest that P. vivax is more virulent than previously thought in tropical endemic areas [5], where genetic diversity is generally high (Table 3). This raises the potential risk of severe malaria in South Korea if genetic diversity further increase to a level observed in tropical endemic areas. Recently, drug-resistant P. vivax has been reported in South Korea [31], which emphasizes increasing importance of assessing genetic diversity and population structure in this area. Since results obtained in this study strongly suggest an increase in the genetic diversity of P. vivax population in South Korea, monitoring of the parasite genetic diversity in this country as well as in neighbouring endemic areas will provide valuable information for developing malaria control strategies in South Korea.