An Evolutionary Arms Race Between Burkholderia pseudomallei and Host Immune System: What Do We Know?

A better understanding of co-evolution between pathogens and hosts holds promise for better prevention and control strategies. This review will explore the interactions between Burkholderia pseudomallei, an environmental and opportunistic pathogen, and the human host immune system. B. pseudomallei causes “Melioidosis,” a rapidly fatal tropical infectious disease predicted to affect 165,000 cases annually worldwide, of which 89,000 are fatal. Genetic heterogeneities were reported in both B. pseudomallei and human host population, some of which may, at least in part, contribute to inter-individual differences in disease susceptibility. Here, we review (i) a multi-host—pathogen characteristic of the interaction; (ii) selection pressures acting on B. pseudomallei and human genomes with the former being driven by bacterial adaptation across ranges of ecological niches while the latter are driven by human encounter of broad ranges of pathogens; (iii) the mechanisms that generate genetic diversity in bacterial and host population particularly in sequences encoding proteins functioning in host—pathogen interaction; (iv) reported genetic and structural variations of proteins or molecules observed in B. pseudomallei—human host interactions and their implications in infection outcomes. Together, these predict bacterial and host evolutionary trajectory which continues to generate genetic diversity in bacterium and operates host immune selection at the molecular level.


INTRODUCTION
Melioidosis is a serious and often fatal neglected tropical infectious disease caused by Burkholderia pseudomallei, an intracellular bacterial pathogen and also ubiquitous in the environment (Holden et al., 2004). Human hosts can acquire the bacterium after direct environmental exposure either through dermal puncture, ingestion of contaminated food or water supplies, or inhalation of contaminated soil or water aerosols. Following an acquisition, B. pseudomallei can replicate in host non-phagocytic and phagocytic cells as well as spreading intracellularly (Willcocks et al., 2016). The bacterium can either kill the human host rapidly (acute melioidosis) or hide within the host body for a long period (chronic melioidosis). With rare exceptions (Abbink et al., 2001;Ralph et al., 2004;Aziz et al., 2020), B. pseudomallei is not known to transmit from person to person; indicating that B. pseudomallei has not evolved virulence mechanisms through consecutive passage in human hosts. Moreover, B. pseudomallei within-host evolution observed in chronic infections are often linked to attenuated virulence (Price et al., 2013;Viberg et al., 2017;Pearson et al., 2020), which is possibly mediated by host immune evasion. While environment exposure can lead to melioidosis, most patients are elderly or have one or more underlying health risks. The most common health condition is diabetes mellitus which presents in up to half of cases (Kronsteiner et al., 2019). Therefore, it is likely that B. pseudomallei has acquired virulence genes for human infection before being exposed to the human hosts, and that these virulence factors are most effective when the host has underlying health issues. It is possible that a successful human infection is mediated by bacterial redundant virulence mechanisms that were acquired and maintained in the environment. In this review, we will leverage recent pieces of genomic studies to better understand the evolution of virulence mechanisms in B. pseudomallei and host defense system.

A MULTI-HOST-PATHOGEN INTERACTION
Interactions between free-living amoebae and environmental bacteria have been proposed to provide an evolutionary training platform for intracellular pathogen such as Legionella pnemophila (Park et al., 2020), Listeria monocytogenes (Schuppler, 2014), and B. pseudomallei (Inglis et al., 2000(Inglis et al., , 2003Noinarin et al., 2016). Many free-living amoebae including the genera Acanthamoeba, Dictyostelium, Naegleria, and Paravahlkampfia share natural habitats with Burkholderia bacteria, including B. pseudomallei (Inglis et al., 2000(Inglis et al., , 2003Noinarin et al., 2016). Amoeba-Burkholderia interactions can range from predator-prey to mutualistic relationships. The fate of these interactions is subjected to the species of amoebae host and Burkholderia (genotypes), stages of host, as well as external factors that modify the outcomes. The social amoebae Dictyostelium discoideum life cycle ranges from unicellular amoebae to multicellular slugs where ten of thousands of single-celled D. discoideum aggregate to form a fruiting body. D. discoideum can directly ingest B. pseudomallei and other Burkholderia bacteria as their food; thereby requiring the bacterium to resist host cell phagocytosis, persist inside a unicellular host, and migrate between host cells through the amoeba cytoskeleton. These processes share many similarities with mammalian host infection and likely prime the bacterium for an intracellular lifestyle. The sentinel cells in D. discoideum, which make up approximately 1% of slug cells, can use antimicrobial defense systems similar to those employed by phagocytes in the human innate immune system (Zhang et al., 2016). The sentinel cells are capable of releasing reactive oxygen species to lyse the soil bacterium. Moreover, they can also release an extracellular trap-the reticulated nets of DNA carrying antimicrobial granules-to kill the invading bacteria. The extracellular trap is an ancient host-defense mechanism common among phagocytic cells across vertebrates and invertebrates, thereby providing a training ground for B. pseudomallei infection in mammalian hosts. External factors such as nutrient availability were shown to determine the fate of D. discoideum-Burkholderia relationship. The association is beneficial to both parties when the nutrient is scarce. Under nutrient limited conditions, Burkholderia-associated D. discoideum produced more spores (Brock et al., 2013;DiSalvo et al., 2015) and had an increased uptake of secondary bacteria that can be used as food (Khojandi et al., 2019) which led to a better survival rate than D. discoideum without Burkholderia association. Under this condition, Burkholderia could be detected inside D. discoideum spores, which enable the bacterium to better disseminate. The mechanism underlining the decision to kill or cooperate is unclear.
The species of Burkholderia and genetic variations within the species have been shown to impact the outcomes of amoeba internalization . It is possible that genetic variations in B. pseudomallei could influence the outcomes of human melioidosis. Heritability scores (Finucane et al., 2015;Hou et al., 2019;Speed et al., 2020) can be used to quantify the proportion of variations in the infection outcomes that can be explained by genetic variations in pathogens and hosts. The technique has been successfully applied to different pathogens (Lees et al., 2017) and human infectious diseases. However, the concept has not been widely adopted for B. pseudomallei -host infection. For B. pseudomallei infection in human hosts with no comorbidity, we estimated that 8% of host mortality could be explained by B. pseudomallei genotypes (h 2 = 0.081, SE = 0.050, p = 0.018; see Supplementary Text for method). This new result highlights a moderate proportion of infection outcomes being explained by the bacterial genetics, and suggests that a substantial proportion of the outcomes could be explained by host genetics.
Interactions between pathogens and Homo sapiens have shaped the evolution of modern humans (Karlsson et al., 2014). Currently, less is known about how melioidosis has shaped the human populations, particularly in melioidosis endemic areas. The host immunity could function as a general defense against all invading pathogen mediated by the innate immune system, or a pathogen-specific defense facilitated through an adaptive immune response. For innate immune response, it is possible that interactions between human hosts and other common pathogens may have driven selection on host defense pathways that affect resistance. These may include human interactions with parasites or bacteria that cause malaria (Shimizu et al., 2000;Malaria Genomic Epidemiology Network, 2008;Timmann et al., 2012;Kariuki et al., 2020), tuberculosis (Mahasirimongkol et al., 2012;Thye et al., 2012), or sepsis (Southeast Asia Infectious Disease Clinical Research Network, 2017;Sweeney et al., 2018) in many tropical and subtropical countries, and cholangiocarcinoma (Zou et al., 2014) in Southeast Asia. Using h 2 , previous works estimate that genetic variations in human host can explain between 32 and 52% of infectious disease susceptibilities, depending on the population studied and the infectious agents (Cooke and Hill, 2001;Jia et al., 2019). Candidate gene approaches have determined host risk-and protective genetic markers for melioidosis infection (West et al., 2012(West et al., , 2013Myers et al., 2014;Chaichana et al., 2017;Dickey et al., 2019), many of which were previously identified as gene targets for other common pathogen infections. To cover melioidosis-specific, as well as broad-pathogen host factors, systemic and reproducible genomewide scans are yet to be conducted.

NATURAL SELECTION IN GENES REQUIRED FOR HOST-PATHOGEN INTERACTION
Given an increasing availability of B. pseudomallei and host genomes, a genome-wide scan for signatures of selection, their directions (purifying, balancing, and positive selection) and their magnitudes could be very useful to investigate host-pathogen protein-protein interactions. As both host and B. pseudomallei migrate, detection of directional selection can be challenging after a population bottleneck. A bottleneck leads to reduced levels of genetic variation and a subsequent loss of a selection signal. For a soil microbe like B. pseudomallei, bacterial migration was shown to be infrequent, with major movements tracked and dated (Pearson et al., 2009;Price et al., 2016;Chewapreecha et al., 2017). The bacterial genome evolution is correlated with a strong geographical signal, with the highest diversity observed in Australia and reduced diversity following dissemination out of the Australian origin (Chewapreecha et al., 2017). Nevertheless, frequent recombination observed in B. pseudomallei can introduce new genetic variations, thereby enabling detection of selection signatures in other contemporary populations outside Australia (Chewapreecha et al., 2019). Natural selection on modern humans as a result of pathogen encounter and migration has been reviewed elsewhere (Karlsson et al., 2014;Sironi et al., 2015;Slodkowicz and Goldman, 2020).
For bacteria, evolutionary pressures on orthologous proteins can be quantified using the ratio between substitution rates at non-synonymous (dN) sites, which could have experienced selection, and synonymous (dS) sites which are presumably neutral. The dN/dS ratio is likely to be more than 1 if pressures favor changes in the protein sequence (positive selection). It is likely to be less than 1 if selections suppress changes in the protein sequence (purifying selection). Several studies (Yu et al., 2006;Losada et al., 2010;Nandi et al., 2010;Hayden et al., 2012;Chewapreecha et al., 2019) and this review utilized dN/dS to estimate selection pressures on coding sequences of the pathogenic B. pseudomallei as well as non-pathogenic but closely related species such as B. thailandensis. Although these studies focused on bacterial populations isolated from different melioidosis endemic regions, performed on different sample size and genetic diversity, they largely highlighted three similar patterns (Figure 1). The first common pattern is an elevated level of positive selection in Burkholderia accessory genes (genes that are variably present across different isolates) compared to core genes (genes that are consistently present in all isolates). This observation highlights the role of accessory genomes in mediating B. pseudomallei adaptation, a feature reflected by the open pangenome and large repertoire of accessory genes observed in B. pseudomallei and B. thailandensis (Nandi et al., 2015;Spring-Pearson et al., 2015;Chewapreecha et al., 2017). The second common pattern is the set of genes under FIGURE 1 | Selection pressure acting on B. pseudomallei population. The histogram summarizes ranges of dN/dS calculated from predicted coding sequences from a collection of diverse B. pseudomallei population from northeast Thailand (Chewapreecha et al., 2019), and B. thailandensis genomes from the public database. B. pseudomallei and B. thailandensis have highly plastic genomes comprising of at least two chromosomes of ∼7-8Mb in size when combined. Using a pan-genome approach, all coding sequences could be categorized as "core" (present in all genomes) or "accessory" (variably present across studied genomes). Accessory genes display an elevated level of dN/dS which is signatures of positive selection or more relaxed purifying selection.
purifying selections in B. pseudomallei and B. thailandensis involved in the bacterial replication, transcription and translation machinery; highlighting the conservation of these genes in both species. The third common pattern is that many genes under positive selection in B. pseudomallei and B. thailandensis are required for environmental survival and exchange of genetic materials including secretion systems, response regulator, heat shock proteins, and integration and/or restriction of horizontal gene transfer. However, genes required for host cell invasions such as cell adhesion, fimbriae, intracellular multiplication and macrophage killing are only positively selected in B. pseudomallei but not B. thailandensis. These distinct signals for positive selection in B. pseudomallei and B. thailandensis could mark a different pathogenic potential of the two species.
Genes under positive selection could primarily offer B. pseudomallei advantages to survive harsh environmental conditions. Given its tropical and sub-tropical habitats, B. pseudomallei experiences drought and heavy rainfall on a seasonal basis. The bacterium can be found in deeper soil layers during the dry season and moves to the soil surface during the rainy season following water movement (Thomas et al., 1979;Currie and Jacups, 2003;Manivanh et al., 2017), highlighting a spatiotemporal transition. The bacterium is also subjected to heat stress. An in vivo study showed that B. pseudomallei had a slower growth rate and a shift in gene expression toward heat-shock proteins and bacterial motility when exposed to temperature stress (Paksanont et al., 2018). Intriguingly, B. pseudomallei is capable of persisting in a nutrient-free distilled water (Pumpuang et al., 2011). The bacterium was subjected to distilled water since 1994 and this experiment is still ongoing, thereby marking an unusual ability of B. pseudomallei to survive in low nutrient media and high osmotic pressure. Lipopolysaccharide, a known virulence factor, was proposed to facilitate B. pseudomallei survival in water. Statistical analysis suggested that there is a correlation between the presence of LPS and rainfall, in particular, LPS serotype B (Shaw et al., 2019). An experimental study further reported that a gene of inner core LPS biosynthesis cluster (waaE: BPSL2510) is vital for a long-term water incubation (Moore et al., 2008). Moreover, an adhesin BPSL1661 was identified as a hub for co-evolutionary signals in B. pseudomallei population (Chewapreecha et al., 2020). The gene was functionally characterized as essential during nutrient starvation (Chewapreecha et al., 2020), thereby highlighting nutrient limitation as a major evolutionary pressure experienced by this microorganism. The presence of B. pseudomallei is strongly associated with the low-nutrient soil (Limmathurotsakul et al., 2010b;Hantrakun et al., 2016). This led to the hypothesis that B. pseudomallei might have a competitive advantage over benign soil microbes in nutrient-depleted soil but is outcompeted in a nutrient-rich environment. The occurrence of B. pseudomallei in the nutrient-poor agricultural field has linked to crop residue burning, a common practice in the tropics where the burning practice could deplete the soil organic matter (Bot and Benites, 2005). An ability to persist in the hostile environmental conditions allows the bacterium to survive before it could seek shelters in a more protecting reservoir inside the single or multicellular hosts.
Genes that facilitate bacterial transition from the environment to hosts are also under positive selection. A genome-wide association study (GWAS) recently identified the bacterial genetic factors that distinguish between B. pseudomallei isolated from the environment and those causing disease in human (Chewapreecha et al., 2019). Genetic variants associated with cell entry and toxin were found to be more prevalent in disease-than environmental isolates, while variants involved in malfunctional cell adherence were found at higher frequency in environmental-than disease isolates. Together, this highlights significant roles of cell adhesion and cell entry in allowing the bacterium to switch to an intracellular lifestyle, either in an amoeba or in a human host where it causes melioidosis. The same study also quantified the numbers of time the variants were gained or lost from the population phylogenetic tree and highlighted multiple gain-and-loss events for both disease-and environmentassociated variants. This genetic evidence suggests a process into which B. pseudomallei can adapt to colonize multiple niches by exploiting existing variations or exchanging disadvantageous for advantageous alleles that promote its survival in a new niche, including the human host.

MECHANISMS THAT GENERATE GENETIC AND MOLECULAR VARIATIONS
Genetic information in both B. pseudomallei and human host can be passed on vertically, although mistakes during replication process could result in small-scale genetic variations. In B. pseudomallei, larger-scale genetic variations can be introduced by recombination event or horizontal gene transfer (HGT). The former process is akin to sexual reproduction and meiotic crossing-overs in its human hosts. Due to contrasting shortand long-generation time of B. pseudomallei (49 min in the log phase; Ou et al., 2005) and host (22-33 years for modern humans; Ewbank, 2016); the host cannot solely rely on new genetic diversity generated when the new offspring is born. A healthy human possesses a large and dynamic repertoire of B cell-and T cell receptors that recognize the invading pathogen and produce a repertoire of antibody that recognizes a variety of antigenic structure. This molecular diversity is achieved without a change in genetic content. The mechanisms that generate the genetic and molecular diversity of the human immune system have been reviewed elsewhere (Nikolich-Zugich et al., 2004;Robinson, 2015;Dendrou et al., 2018;van den Broek et al., 2018;Adams et al., 2020).
Small-scale genetic variations in B. pseudomallei can be introduced by point mutation which substitutes one nucleotide with another or microindel (an insertion or deletion) that impact 1-50 bp. B. pseudomallei substitution rate was reported to be 1.7-4.9 × 10 −7 substitutions per site per year (Viberg et al., 2017;Pearson et al., 2020), which is comparable to those detected in other Burkholderia genera (Lieberman et al., 2014). The rates of indels have not been quantified in B. pseudomallei. However, indels can be detected from isolates obtained from acute and chronic host infections (Hayden et al., 2012;Limmathurotsakul et al., 2014a;Viberg et al., 2017) as well as from the environment (Rachlin et al., 2020), thereby highlighting the role of indels in B. pseudomallei evolution.
Medium-scale genetic variations can be brought in by homologous recombination, with each event contribute to a median recombining size of 5 kb (Nandi et al., 2015) (range 3 bp to 71 kb). A single recombination event can introduce 7.2 times greater nucleotide polymorphisms than a single substitution event (average r/m = 7.2). However, the amounts of SNPs being introduced per recombination event is subjected to the genetic distance between DNA donor and recipient. Sources of imported DNA could be from the same species donor, closely related species such as B. thailandensis and/or other soil microbes. Different B. pseudomallei lineages have been shown to recombine at different rates (Nandi et al., 2015), a characteristic also observed in other recombinogenic bacteria (Chewapreecha et al., 2014;Sanchez-Buso et al., 2014;David et al., 2017). Interestingly, not all DNA donorrecipient pairs are possible, suggesting a structure to the genetic flux within B. pseudomallei population. A lineagespecific restriction-modification system was shown to act as a barrier that restricts gene flow (Nandi et al., 2015). These modification control systems are based on DNA methylation which allows the bacterium to discriminate between correctly methylated "self " DNA, and inappropriately methylated or unmethylated "non-self " DNA. Correctly methylated DNA can be taken up and integrated into the bacterial genome, whereas inappropriately methylated or unmethylated DNA will be degraded. Interestingly, many recombination "hotspots" are focused on bacterial virulence factors (Nandi et al., 2015), highlighting the role of recombination in tuning the bacterial ability to infect hosts.
Large-scale genetic variations can be introduced by HGT or large-scale insertion or deletion. HGT includes the uptake of foreign DNA that can subsequently be integrated into the chromosomes and resulted in the regions termed genomic islands (GIs) (Koonpaew et al., 2000;Holden et al., 2004;Duangsonk et al., 2006). The repetitive genetic composition of B. pseudomallei genomes such as tandem repeats and transposons is believed to aid the integration of horizontally acquired elements into the chromosomes (Holden et al., 2004;U'Ren et al., 2007). Despite several plasmid elements being identified in B. pseudomallei genomes (Holden et al., 2004), four complete plasmids has been characterized so far (NCBI, December 2020) (Nandi et al., 2015). With long-read sequencing technology, the plasmid could be more reliably assembled, thereby allowing more genetic variations to be studied. At least 16 GIs had been identified in the K96243 genome (Holden et al., 2004). Several studies highlight distinct geographical distribution of GIs or their combinations (Ou et al., 2005;Duangsonk et al., 2006;Sim et al., 2008;Tuanyok et al., 2008), leading to a hypothesis that each GI or combination of GIs may confer different fitness under different environments (Sim et al., 2008). GIs carry several virulence factors, many of which are known to interact with hosts, including filamentous hemagglutinin adhesin (fhaB3), and a member of two-partner secretion system (bpaAB) (Sim et al., 2008).

GENETIC, STRUCTURAL OR MOLECULAR VARIATIONS OBSERVED IN HOST-PATHOGEN INTERACTION
Genetic, structural and molecular variations and their functional impact to B. pseudomallei -host interaction have been the center of melioidosis research over the past few decades.
We will summarize genetic and structural diversity in bacterial genes known to mediate the infection (Figure 2 and Supplementary Table 1). The genes discussed here are shown to be essential for infection by a genome-wide saturation mutagenesis (Moule et al., 2015) and also expressed during in vivo infection (Ooi et al., 2013). When host partners are known; variations detected in host proteins, and the host response are also described (Figure 3). It should be noted that there are large overlaps in the host immune system and the examples described here only represent a fraction of the whole machinery.
Genetic Variations in B. pseudomallei Flagella Systems and Host Toll-Like Receptor 5 (TLR5) Flagella systems of B. pseudomallei enable the motility of the bacterium intracellularly as well as in the environment (Figure 2; DeShazer et al., 1997;French et al., 2011). Flagella filaments are constituted from thousands of flagellin protein monomers, encoded by fliC gene, strung into protofilaments before being braided to form a flagellum (Samatey et al., 2001). Flagellin has been known for its extreme diversity with 113,285 unique nucleotide sequences across the prokaryote phyla (Hu and Reeves, 2020). Flagellin inactivation in B. thailandensis resulted in the reduction of bacterium intracellular motility and cell-to-cell spread (French et al., 2011). Using PCRrestriction fragment length polymorphism analysis, four different alleles of flagellin protein (fliC: BPSL3319) were identified from 100 Malaysian clinical isolates (Tay et al., 2010), and more is expected if investigated with nextgeneration sequencing. Proteins in the flagella systems were reported to elicit host immune response, highlighting them as melioidosis vaccine candidates (Brett et al., 1994;Chua et al., 2003;Chuaygud et al., 2008;Gregory et al., 2015;Koosakulnirand et al., 2018).
In addition to flagella, a fimbrial gene cluster which is required for cell adherence also displays genetic variation with distinct alleles being predominantly detected in different geographical locations. A yersinia-like fimbrial (YLF) gene cluster, believed to be horizontally acquired was shown to be more prevalence in Southeast Asia and thus was used as a marker for the introduction of B. pseudomallei from this region (Tuanyok et al., 2007;Sarovich et al., 2014;Chewapreecha et al., 2017). A putative type I fimbrial protein BPSL1626 was shown to elicit an immune response and have a potential as a vaccine candidate against melioidosis (Capelli et al., 2018). Whereas, B. thailandensis-like flagellum and chemotaxis (BTFC) gene cluster is believed to be an ancestral sequence in B. pseudomallei and is most common in Australia (Tuanyok et al., 2007). The YLF and BTFC gene clusters are mutually exclusive between the two endemic areas where the latter may implicate in the spread between cell-to-cell by the flagellar protein leading to the formation of multinucleated giant cells (MNGCs) and eventually apoptosis/cell death (French et al., 2011).  (Ooi et al., 2013), with a subset that reported genetic variations marked in red. The shape and location of each individual gene indicate the gene function and cellular compartment, respectively. All annotated genes or operons and their functions are described. The bacterium displays a repertoire of antigenic variations, including lipopolysaccharides (LPS), capsular polysaccharides (CPS) and surface proteins. B. pseudomallei LPS is immunologically classified into a number of serotypes A, B, and B2; with each serotype reported to be heterogeneously distributed across distinct geographical locations. Another highly diverse virulent protein is a fimbrial protein which displays a strong geographical distribution between Australia and Asia. Strains from Asia commonly possess a Yersinia-like fimbrial (YLF) gene cluster that believed to be horizontally acquired. B. pseudomallei carries 4 different types of CPS: CPS I, CPS II, CPS III, and CPS IV. A full gene description is provided in Supplementary Table 1. FIGURE 3 | A summary of human host immune components used in defending against B. pseudomallei. Genes or molecules that show genetic variations are marked in red. Upon infection, the bacterial antigens are recognized by host receptors including Toll-like receptor (TLR) and HLA which also display large genetic diversity. The latter is reported to be varied by ethnic groups (de Bakker et al., 2006;Gourraud et al., 2014). Following an entry into host, B. pseudomallei faces the innate axis of the immune system. Macrophages and neutrophils are recruited to the site early upon the infection; these cells are reported to be essential to the early bacterial containment and clearance Easton et al., 2007), though excessive recruitment of neutrophils may have a negative outcome that allows B. pseudomallei to propagate intracellularly (Ceballos-Olvera et al., 2011). These innate phagocytes possess pattern recognition receptors (PRRs), such as surface receptors TLR2, TLR4, and TLR5, as well as, cytosolic receptors NOD2; these are reported to be vital for the fight against B. pseudomallei (West et al., 2008;Myers et al., 2014;Weehuizen et al., 2015;Birnie et al., 2019). The signals transduced by these receptors result in mobilization of nuclear factor NF-κB which trigger appropriate immune responses including the synthesis of pro-inflammatory cytokines and initiation of the downstream adaptive immune cascades (Pothlichet and Quintana-Murci, 2013). It is through the antigen-presenting cells and their appropriate antigen-HLA class II complexes that enable the activation of CD4 + T lymphocytes. The naïve CD4 + T cells sit at the central part of the adaptive axis. They can differentiate into Th1 cells which facilitate cell-mediated immune response by CD8 + cytotoxic lymphocytes. A study has shown that strong CD4 + and CD8 + T cell response was elicited during acute melioidosis and the lower cellular response was correlated to fatality (Jenjaroen et al., 2015). On the other hand, naïve CD4 + T cells can also mature to be Th2 cells which initiate class switching of B cells and support humoral immune response.
It has been established that flagellin is a ligand of TLR5 (Figure 3). The recognition stimulates pro-inflammatory responses including the rise in intracellular calcium ion and upregulation of pro-inflammatory cytokine TNF-α and IL-6 (Hayashi et al., 2001;Chen et al., 2007). TLR5 stop codon polymorphism, TLR5 1174C>T or rs5744168 was strongly associated with protection against fatality as well as organ failure in a case-control cohort study. When challenged with B. pseudomallei, TLR5 1174C could mediate the activation of NF-κB, while TLR5 1174T could not; additionally, TLR5 1174T saw reduced flagellin-induced cytokines levels (West et al., 2013). Another independent investigation on the same population also found the association between the truncated TLR5 1174C>T variant and survival from acute melioidosis, as well as, a lower rate of bacteraemia (Chaichana et al., 2017). Furthermore, both studies also found TLR5 1174C>T variant with a lower level of anti-inflammatory IL-10, which the authors suggested the possibility in which the mortality risk may be modulated by TLR5-driven IL-10 release. Interestingly, a recent study has suggested that the effect of this hypofunctional TLR5 variant may not be restricted to flagellin-driven pathway (Dickey et al., 2019). Another polymorphic variant is TLR5 1846T >C which was also associated with protection against death and blunted flagellin-driven cytokine response; however, the authors also reported high linkage disequilibrium of the variant with TLR51174C>T , which might reduce the confidence of causal relationship.

Structural Variations in B. pseudomallei LPS and Genetic Variation in Host Toll-Like Receptor 2 and 4 (TLR2 and TLR4)
B. pseudomallei possesses an extensive network of polysaccharides on its outer membrane, namely capsular polysaccharide (CPS) and lipopolysaccharide (LPS) (Figure 2, for a full review on CPS and LPS in Burkholderia spp.; Cloutier et al., 2018). Both are known to play a vital role in virulence of melioidosis and have been used as subunits for vaccine development (Brett and Woods, 1996;Nelson et al., 2004;Wikraiphat et al., 2009;Scott et al., 2014b). Based on a cellular compartment, B. pseudomallei LPS can be divided into (Knirel et al., 1992): lipid A-an endotoxic component embedded in the phospholipid bilayer of the outer membrane; inner and outer core oligopolysaccharide; and O-antigen (Novem et al., 2009). Multiple in vitro studies showed LPS challenge could mount innate and adaptive immune responses, and Nitric Oxide (NO) production, with different LPS serotypes (A, B, and B2) reported to mount different magnitude of responses (Norris et al., 2017). Some genetic data of LPS biosynthesis operons were available. However, many of these were generated from short-read sequencing platforms, which are not ideal to investigate LPS due to the repetitive nature of the LPS locus. Long-read technologies, on the other hand, could overcome this assembly issue, giving complete contigs of the bacterium genome.
To date, 4 serotypes of LPS were previously characterized. Type A LPS is a majority serotype found in both endemic areas, namely Southeast Asia and Australia (Anuntagool et al., 2000(Anuntagool et al., , 2006Tuanyok et al., 2012). Whereas, serotype B is less common in the endemic areas (Australia and Southeast Asia) and is prevalent in both clinical and environmental origins (Tuanyok et al., 2012). However, a recent analysis of clinical B. pseudomallei isolates revealed that Serotype B was highly predominant in India (Shaw et al., 2019). LPS serotype B2 is a variant of LPS serotype B and more commonly found in Australia and Papua New Guinea (Tuanyok et al., 2012). In contrast, LPS R or Rough serotype lacks O-antigen moiety of the LPS structure. It was identified using SDS-PAGE and silver staining techniques with no O-antigen ladder pattern. Type R is relatively rare and found only from the Australian clinical and environmental strains (Anuntagool et al., 2006). Interestingly, it is frequently prevalent in the patients with a relapse history of melioidosis, however, there is no direct association between them (Limmathurotsakul et al., 2014b).
A comparative genomic analysis revealed that there are distinct variations in the core compositions of the O-antigen LPS biosynthesis gene clusters between A, B, and B2 serotypes. However, there are several genes conserved among them. In addition, type A LPS gene operon is also observed in both B. thailandensis and B. mallei, closely related Burkholderia species (Tuanyok et al., 2012). Interestingly, genomic analysis of clinical isolates from Madagascar revealed that there is a 13.5 kb deletion observed in the LPS biosynthesis gene cluster of serotype B, conserving only some genes in the cluster that are essential for the biosynthesis of LPS B2 . The lack of some core genes in the LPS biosynthesis cluster could lead to the reduction of its serological properties. This is supported as some evidence suggested that the strains with LPS type B2 become sensitive to 30% normal human serum whereas the strains with LPS type B remain resistant (Tuanyok et al., 2012).
On the other hand, a single nucleotide insertion of the wbiI gene is observed in the LPS biosynthesis gene cluster of a rough serotype in a patient with more than 16-year chronic lung infection associated with melioidosis. This frame-shift mutation of wbiI gene (an essential gene for the O-antigen synthesis) disrupts the epimerase/dehydratase function of this gene and results in the loss of O-antigen moiety, possibly switching the serotype of B. pseudomallei isolates from type A to type R (Price et al., 2013;Pearson et al., 2020). Pearson and colleagues also revealed the nucleotide insertion of D512_15771 (wbiH) and D512_06755 gene are observed in B. pseudomallei MSHR6686 of the same patient. The mutation in both genes may confer to the reduction in LPS modification and production, assisting in the escape from the host immune response. In addition, a partial deletion of D512_20407, wbiA homolog is also identified in MSHR6686 (Price et al., 2013;Pearson et al., 2020). However, genomic analysis of the rough serotype that is naturally found in the environment or initial infection has not been done to determine the genetic makeup of this serotype.
Not only genetic heterogeneity exists among the serotypes of LPS, further structural diversity is also observed at the O-antigen of LPS serotypes Norris et al., 2017). The O-antigen is one of the LPS components with a structure of unbranched disaccharide repeat units. Remarkably, the structural modifications of these sugar chains were observed where substitutions of 2-O-methylated and 4-O-acetylated at talose residues were observed only in about 33% of the LPS serotype A. Whereas the rest bear 2-acetyl substituents at the same residues (Perry et al., 1995). More recently, a structural analysis of O-antigen in the serotype A reported that the modification of the talose residue is more complex than what was previously reported (Heiss et al., 2013). Although multiple gene inactivation studies revealed that wbiA and oacA genes are essential for these modifications at the talose residues (Brett et al., 2003, it is possible that structural diversities of the LPS O-antigen are further modified at post-transcriptional and posttranslational levels, for instance, a length variation of O-antigen chain observed in Escherichia coli O9a which is tightly controlled by the biosynthetic enzymes dynamic (King et al., 2014).
The LPS is common across Gram-negative bacteria and also a well-established pathogen-associated molecular pattern (PAMP) that can trigger the pro-inflammatory innate immune response, such as translocation of NF-κB and TNF-α cytokine release, via its interaction with TLR4, CD-14 and adaptor protein MD-2 (Park and Lee, 2013). Although evidence has been inconclusive, LPS sensing was shown to be through TLR4 in murine models, and TLR4 as well as TLR2 (Wiersinga et al., 2007;West et al., 2008) in human models (Figure 3; Weehuizen et al., 2015). Polymorphisms of TLR4 were reported in humans with TLR4 1196C>T allelic variant associated with protection against melioidosis when compared to non-hospitalized controls (West et al., 2012). In addition, TLR4 rs10828066 SNP variant was significantly associated to be protective against melioidosis in both adjusted non-hospitalized and hospitalized control groups, whereas TLR4 rs960312 was associated with bacteraemic or lung melioidosis (West et al., 2012).

Genetic Variations in B. pseudomallei CPS and Host B Cell Repertoire
CPS is one of the key virulence determinants in B. pseudomallei (Figure 2). To date, four different types of CPS (CPS I-IV) have been described (Holden et al., 2004). CPS-I biosynthesis cluster is a large 34.5 kb operon harboring on the chromosome one of B. pseudomallei (wcbT to manC). A study in an animal model reported that the CPS-I cluster is required for the full virulence of B. pseudomallei (Sarkar-Tyson et al., 2007). CPS-I exhibits an immunogenic role and the absence of this gene cluster resulted in the attenuation of B. pseudomallei in a mouse model (Atkins et al., 2002;Reckseidler-Zenteno et al., 2005;Parthasarathy et al., 2006). A passive protection was observed when murine intranasal infection models were immunized with anti-CPS monoclonal antibodies (AuCoin et al., 2012), with the later study reported that conjugated CPS could provide the highest degree of protection (Scott et al., 2014a).
The structural basis of CPS-I is highly conserved (Perry et al., 1995). Remarkably, CPS-I is believed to be horizontally transferred between species as the GC content of this CPS-I cluster is about 58% (Reckseidler et al., 2001). In addition, it is possible that the horizontal acquisition of this gene cluster was a key event in the pathogenic evolution of B. pseudomallei compared to B. thailandensis where the cluster is absent (Yu et al., 2006). Nevertheless, several genes involved in the biosynthesis of CPS-I of B. pseudomallei were reported in B. thailandensis with a relatively low sequence identity (75%). In addition, a novel variant of B. thailandensis strains (BTCV) is found to carry a B. pseudomallei like CPS-I (95% sequence identity) gene cluster with an identical organization. Although the origin of the CPS-I gene cluster in B. pseudomallei has yet been identified, several genes involved in the sugar biosynthesis appear to be orthologs to the genes identified in the Yersinia species (Y. pseudotuberculosis H892/87) and some in other gram-negative bacteria (Cuccui et al., 2012). It is possible that the CPSI gene cluster is not entirely acquired from one bacterial species but several organisms. Population studies also detected genetic variations in CPS-I (Challacombe et al., 2014;Chewapreecha et al., 2017). Whether this genetic variation has implications on virulence is a subject of further investigation. Another CPS-I variant could be found in B. pseudomallei isolated from a 16-year chronic melioidosis patient. A single nucleotide insertion is observed within wcbR, an important component involved in the fatty acid synthesis of CPS-I, causing a frameshift mutation in the earlier samples isolated from this patient (Price et al., 2013;Pearson et al., 2020). A CPS-I deletion region, which includes wcbR gene, has shown in the reduction of CPS production in hamster model but not absence entirely (Gutierrez and Warawa, 2016). Price et al. (2013), suggested that this degree of CPS-I-dependent virulence decreases and may consequently be a critical step in the progression for melioidosis to become a chronic-carriage disease.
Unlike CPS-I, CPS-II and CPS-III have been associated with B. pseudomallei persistence in the environment (Reckseidler-Zenteno et al., 2009, while CPS-IV is less wellcharacterized. Up to date, no reports suggest the genetic variations of CPS-II, -III, and -IV biosynthesis clusters in B. pseudomallei. The diversity and function of these capsules in the pathogenesis, immunomodulation and environmental adaptation of B. pseudomallei warrant further studies. We next considered host partners that interact with B. pseudomallei CPS. The CPS has been shown to elicit strong host immune response (Reckseidler et al., 2001;Atkins et al., 2002). CPSs from various B. pseudomallei strains were recognized by the same group of monoclonal antibodies which suggests a limited number of epitopes of this molecule (Zou et al., 2008). The molecule's monomers can cross-link to B cell receptors and, at an appropriate density, induce a downstream humoral response that is distinct from T cell-dependent pathway (Clarke et al., 2013;Akkaya et al., 2020). Antigen recognition through B-cell receptors is formed through random somatic changes of germline DNA. This results in a repertoire of distinct sequences that enable antigen recognition across wide ranges of pathogens. The B-cell receptors dynamic in melioidosis has not been studied.

Genetic Variations in B. pseudomallei Adhesins and Autotransporters
Filamentous Hemagglutinin Adhesin (FHA) of B. pseudomallei is highly diverse (Figure 2). Early investigation of genomic islands (GIs) from five B. pseudomallei strains identified three different fhaB gene clusters on different GIs of the bacteria. The bacterium could carry multiple fhaB gene clusters, many of with carrying either a combination of cluster I (GI5a/GI5a.1/GI5a.2) and cluster III (GI16/GI16.1), or cluster III alone . More genetic variations of the FHA loci were reported from the Australian isolates (Chewapreecha et al., 2017). The fhaB3 gene (BPSS2053) has been characterized as an important virulence factor of B. pseudomallei; enabling bacterial binding to the host epithelial cells (Sim et al., 2008), and an antimacrophage factor (Dowling et al., 2010). It was found in all isolates from Thailand but found in only 83% of Australian strains where the absence of this gene in Australia population correlated with the skin abscess formation and lower mortality rate (Sarovich et al., 2014).
Autotransporter (AT) proteins are one of the largest family of the secretion systems in Gram-negative bacteria, allowing the transportation across the bacterial membrane as well as involving the virulence and immunogenicity of the pathogens (for review Lazar Adler et al., 2011). In B. pseudomallei K96243, a sequence analysis revealed to have at least 11 ATs located in the genome, including putative Trimeric Autotransporter Adhesins (TAAs), bimA (Stevens J.M. et al., 2005) and boaB, twopartner secretion system (TPS), bpaA, bpaB, and bpaD (Campos et al., 2013). Inactivation in some of these genes attenuated and reduced the intracellular survival of B. pseudomallei in the macrophage-like cells (Lazar Adler et al., 2015). ATs are variably present in the genomes. The gene boaB (BPSL1705), which aids bacterial binding to host respiratory cells (Balder et al., 2010), was found to be absent in several strains in Africa, Brazil . When ATs are present, genetic variants can also be observed. Australian-specific and Africanspecific variants have been reported in bpaA, bpaB, and bpaD (BPSS1434, BPSL2063, and BPSS0088, respectively) (Brown et al., 2004;Tuanyok et al., 2008).
Several experiments demonstrated the role of BimA in actin polymerization and motility (Stevens M.P. et al., 2005;French et al., 2011), promoting a movement of cell-to-cell spread and within the host cells (Lazar Adler et al., 2011). Furthermore, BimA is antigenic (Suwannasaen et al., 2011) and seroactive (Felgner et al., 2009). Two variants of bimA (BPSS1492) were identified: B. pseudomallei bimA (bimA BP ), and bimA-like B. mallei (bimA BM ). The latter is an ortholog of bimA from Burkholderia mallei (95% sequence identity) and displays the same domain organization: a single actin monomer binding motif (WH2), a proline-rich domain and a transmembrane anchor domain (Sitthidet et al., 2008). Whereas the bimA BP consists of two predicted WH2 domains, a proline-rich domain, a membrane anchor domain and an additional predicted casein kinase II (Stevens M.P. et al., 2005). Of note, the variation was also observed between the predicted proline-rich domain of the bimA BM and B. mallei bimA where the former has fewer motifs (Sitthidet et al., 2008). Distinct geographical distribution of bimA BP and bimA BM have been reported (Sitthidet et al., 2008;Sarovich et al., 2014;Shaw et al., 2019). Benanti et al. (2015) illustrated that both BimAs are functionally similar and associated with the nucleation and elongation of actin filaments with more plaque formation observed in bimA BP variant compared to bimABM. A subunit vaccine investigation reported a significant increased survival rate of BALB/c mice immunized with BimA BM recombinant, when the mice were challenged with both B. mallei and B. pseudomallei (Whitlock et al., 2010).

Genetic Variations in B. pseudomallei Secretion Systems
Several secretion systems and secreted proteins in B. pseudomallei are vital to combat the host immune defense and environmental stresses. B. pseudomallei has at least three type III secretion systems (T3SSs) (for review, Vander Broek and Stevens, 2017)and six type VI secretion systems (T6SSs).
Among T3SSs, T3SS-3 is better characterized and believed to be an integral part of the full B. pseudomallei virulence in mice and hamster models, as well as associated with the intracellular survival and dissemination. K-mer based sequence analysis has recently demonstrated that variations in the T3SS-3 gene cluster were detected across the global population of B. pseudomallei (Chewapreecha et al., 2017). This included bsaU(BPSS1539), bsaR (BPSS1542), bsaP (BPSS1544), and bsaN (BPSS1546) (Chewapreecha et al., 2017) which are believed to be involved in intracellular escape (Pilatz et al., 2006), predicted chaperone protein (Panina et al., 2005), T3SS-3 secretion regulator (Broek et al., 2015), and T3SS-3 regulator (Chen et al., 2014), respectively. In addition, a partial deletion of genes in the T3SS is also observed in a patient with persisted infection of melioidosis as a result of within-host adaptations (Pearson et al., 2020). For T6SSs, T6SS-5 is better characterized (for review, Lennings et al., 2019)and functionally confirmed to mediate the translocation of effector proteins via contact-dependent manner (Silverman et al., 2012) and MNGC formation Chen et al., 2011). In addition, an experimental study suggested that T6SS-5 may involve in the intracellular survival of B. pseudomallei in macrophages where the expression level of T6SS-5 is dependent on virAG and bprC regulatory gene (Chen et al., 2011). Genetic variations are also observed in several genes clustered T6SS-5 with distinct geographical distribution (Chewapreecha et al., 2017). Remarkably, a large deletion region on chromosome 2 of an environmental B. pseudomallei A4 and an isogenic strain of K96243 has been demonstrated of which they failed to form a plaque in epithelium cells. Of those, the deleted genes are including genes involved in T3SS-3 and T6SS-5 systems. Interestingly, Saiprom et al. (2020) also identified an absence of multiple T3SS-1 genes of previously described Thai environmental strain RF80.

Genetic Variations in Human Leukocyte Antigen (HLA)
For the rest of the bacterial factors that do not get recognized by any specific innate pattern recognition receptors, these antigens will be processed and presented on HLA of antigenpresenting cells such as dendritic cells before getting recognized by appropriate T cell receptors (TCR) of the adaptive immune axis; this HLA-peptide-TCR interaction kicks start the adaptive immune response. There are three different classes of HLA: class I interacts with TCRs on CD8 + T cells while class II binds to TCRs on CD4 + T cells, and the less well-established class III which is not involved in antigen processing and presentation (Dendrou et al., 2018). Acute melioidosis patients with diabetes mellitus were reported to have lower HLA-DR expression on plasmacytoid dendritic cells than the non-diabetic diseased group (Kronsteiner et al., 2019). In addition, in non-diabetic patients, fatal cases presented with significantly lower expression of HLA-DR on monocytes and plasmacytoid dendritic cells, compared to the survived cases (Kronsteiner et al., 2019). Similar observation was found in murine infection models where B. pseudomalleiinfected BALB/c and C57BL/6 mice showed reduced expression of MHC class II on plasmacytoid dendritic cells, however, this was not statistically significant (Williams et al., 2015).
Genetic variations in HLA have been linked to the melioidosis outcomes. A work conducted in an endemic area of Thailand has compared HLA allele frequencies in melioidosis cases and healthy controls. The authors reported a significant increase of DRB1 * 1602 frequency in melioidosis patients, compared to healthy controls. Moreover, an increase in HLA-DRB1 * 1602 and a decrease in HLA-DQA1 * 03 allele frequencies were associated with septicaemic cases of melioidosis (Dharakul et al., 1998). Another study conducted on the same Thai population screened a panel of various HLA class I genotype frequencies in survived and fatal cases of melioidosis. The authors found that HLA-B * 46 and HLA-C * 01 were associated with increased mortality; they were also reported to be in linkage disequilibrium (Dunachie et al., 2017). Interestingly, HLA has been linked with diabetes mellitus, the prominent comorbidity of melioidosis with 12fold increased risk of developing the disease (Limmathurotsakul et al., 2010a). Several recent publications using GWAS have identified variations of HLA and their corresponding protective or predisposing association with type 2 diabetes (Williams et al., 2011;Scott et al., 2017;Zhao et al., 2017). These prompt further characterization of HLA variants and their relationship with melioidosis progression when modulated by patient diabetic status. However, genetic studies on HLA and melioidosis have suffered from small sample size and availability of reliable HLA typing platform. This is largely impeded by the underrepresentation of genomic data from the population from melioidosis endemic areas. The availability of data is crucial as this improves imputations and discovery of new causal variants and disease association.

DISCUSSION
The advancement in omic technologies has improved our understanding of co-evolution between B. pseudomallei and different hosts, thereby guiding better control policy, treatment option and vaccine design. Proteins or molecules that participate in host-B. pseudomallei interaction are extremely variable. Their variations have been seen at the genomic, epigenetic, transcriptomic and proteomic levels. Although genome data for human population from melioidois endemic areas is still scarce, genome data for B.pseudomallei has been accumulating (Holden et al., 2004;Hayden et al., 2012;Price et al., 2013Price et al., , 2016Sahl et al., 2013Sahl et al., , 2016Daligault et al., 2014;Bugrysheva et al., 2015;Chen et al., 2015;Hsueh et al., 2015;Johnson et al., 2015a,b;McRobb et al., 2015;Nandi et al., 2015;Sidjabat et al., 2015;Song et al., 2015;Spring-Pearson et al., 2015;Viberg et al., 2015;Chapple et al., 2016;Aziz et al., 2017;Chewapreecha et al., 2017Chewapreecha et al., , 2019Podnecky et al., 2017;Viberg et al., 2017;Webb et al., 2019). When combined with spatial and temporal information, this allows further exploration of allelic variants and a shift in allele frequency over space and time. Moreover, a genome-wide saturation mutagenesis can aid prediction of essential genes required under certain conditions (Moule et al., 2015). This can be coupled with transcriptome information to understand variations in the expression patterns (Chieng et al., 2012;Ooi et al., 2013;Price et al., 2018) during the course of infection, and across different host types. A dual hostpathogen transcriptome study has not been conducted for melioidosis but is promising to provide valuable insight into the interaction as well as variations that lead to different infection outcomes.
In this article, we mainly explored variations of B. pseudomallei genes implicated in human host infections at the genetic level. Many of these genes display a strong geographical signal which could either be a result of a founder effect following a migration out of Australia, or an acquisition of new alleles required for local adaptation after an introduction to new geographical location. For each virulence gene, we also noted co-existence of multiple alleles in B. pseudomallei population isolated from the same geographical region. The genetic polymorphism could be maintained by balancing selection where each co-existing allele must be favored under different condition. In the context of virulent genes, these could involve competitions with different soil organisms. Signals for positive selection could be detected in genes that promote bacterial survival under hostile environment, and genes required for cell entry and adaptation to an intracellular lifestyle. The latter can be grouped as virulence genes although they may primarily be used in amoeba hosts rather than mammalian hosts. For multi-host-pathogen interaction, it is thus essential to consider virulence in a broad host and environmental context. The picture is far from complete at the moment, but more incoming omic data is promising to shed light on this complex relationship.

AUTHOR CONTRIBUTIONS
CCe conceived the study and performed data analyses. CCo and PB performed literature review. AS performed data analyses. All authors wrote the manuscript.

SUPPLEMENTARY MATERIAL
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb. 2020.612568/full#supplementary-material Supplementary Text | We estimated the proportion of infection outcome variations (survival vs. mortal outcomes after 28 days) explained by B. pseudomallei nucleotide polymorphisms (SNPs). SNPs were called from Chewapreecha et al. (2019Chewapreecha et al. ( , 2020 using kSNP v3 software. A filtering for minimum allele frequency at 0.01 (MAF = 0.01) resulted in 65,379 SNPs for the analysis. B. pseudomallei population was defined by PopPUNK (Lees et al., 2019). We employed the Genome-wide Complex Trait Analysis (GCTA) (Yang et al., 2011) to calculate h 2 using genomic restricted maximum likelihood (GREML) while adjusting for population structure using population clusters defined by PopPUNK.