Immune targets for schistosomiasis control identified by a genome-wide association study of East African snail vectors

Abstract

Schistosomiasis, afflicting >260 million people worldwide, could be controlled by preventing infection of freshwater snail vectors. Intestinal schistosomiasis, caused by Schistosoma mansoni, occurs predominantly in Sub-Saharan Africa and is vectored by Biomphalaria sudanica and related Biomphalaria species. Despite their importance in transmission, very little genomic work has been initiated in African snails, thus hindering development of novel control strategies. To identify genetic factors influencing snail resistance to schistosomes, we performed a pooled genome-wide association study (pooled-GWAS) on the offspring of B. sudanica collected from a persistent hotspot of schistosomiasis in Lake Victoria, Kenya, and exposed to sympatric S. mansoni. Results of the pooled-GWAS were used to develop an amplicon panel to validate candidate loci by genotyping individual snails. This validation revealed two previously uncharacterized, evolutionarily dynamic regions, SudRes1 and SudRes2, that were significantly associated with resistance. SudRes1 includes receptor-like protein tyrosine phosphatases and SudRes2 includes a class of leucine-rich repeat-containing G-protein coupled receptors, both comprising diverse extracellular binding domains, suggesting roles in pathogen recognition. No loci previously tied to schistosome resistance in other snail species showed any association with compatibility suggesting that loci involved in the resistance of African vectors differ from those of neotropical vectors. Beyond these two loci, snail ancestry was strongly correlated with schistosome compatibility, indicating the importance of population structure on transmission dynamics and infection risk. These results provide the first detail of the innate immune system of the major schistosome vector, B. sudanica, informing future studies aimed at predicting and manipulating vector competence.

Keywords: Biomphalaria; GWAS; Schistosoma; amplicon panel; pathogen resistance.

Fig. 1.

Results of the pooled genome-wide association study (pooled-GWAS) identifying variants associated with resistance to Schistosoma mansoni in the Biomphalaria sudanica genome. Fisher’s exact test p values for all pooled-GWAS variants, are arranged horizontally based on contig orthology to 18 chromosomes (x-axis labels) of the B. glabrata genome (xgBioGlab47.1, NCBI RefSeq: GCF_947242115.1) and linkage map analysis (Dataset S2, Fig. S2–S4). Pooled-GWAS dual-variants, defined as two or more proximate (<50 kb and >1.5 kb apart) variants strongly associated with S. mansoni resistance (p ≤ 2.5e-9) are red. All others are blue. All dual-variants and singleton-variants with p ≤ 1e-17 are labeled red and blue, respectively, with their contig and contig position. Unassigned contigs could not be unambiguously mapped and are mostly small and/or repetitive.

Fig. 2.

Amplicon panel validation of pooled-GWAS results. (A) Proportion of ancestry from Biomphalaria sudanica Population 1 (red) or Population 2 (blue) as predicted by ADMIXTURE (31) is correlated with Schistosoma mansoni infection (p < 1e-14), such that majority-Population 1 snails (Group A) are 68% positive while the majority-Population 2 snails (Group B) are 41% positive. Total number of samples included in analysis was 503 including 496 genotypedvalidation and genotyped-pooled-GWAS B. sudanica from this study, 5 inbred B. sudanica previously sequenced (27) and 2 outbred B. choanomphala. (B) Ordered Fisher’s exact test p values per variant of genotyped-validation samples within ancestry groups. SudRes1 variants c6:3,490, c6:31,057, and c6844:72,816 are significant validated-variants (p < 1.3e-04, Bonferroni adjusted significance threshold shown by dotted line) within ancestry Group B (blue dots) and most other top-outliers (1.3e-04< p < 1.0e-03) are in SudRes1 or SudRes2. (C) Ordered additive regression p values per variant of genotyped-validation samples, after accounting for ancestry. SudRes2 variant sc94:2,166,296 is a significant validated-variant (p < 2.3e-04, Bonferroni adjusted significance threshold shown by dotted line) and most other top-outliers (1.3e-04< p < 1.0e-03) are in SudRes1 or SudRes2. (D) Allele and genotype counts for the most significant validated-variants from SudRes1 (c6:3,490) and SudRes2 (sc94:2,166,296), which were both dominant markers, and the representative codominant marker variants from SudRes1 (c582:65,596) and SudRes2 (sc94:2,174,117), for S. mansoni infection positive (+) and negative (−) snails. Values are read counts for each allele in the pooled-GWAS (“GWAS”), and genotype counts in genotyped-validation snails (“Validation”). (E) Means and standard errors for each variable in the best-fitting multiple regression model. Ancestry is the proportion of Population 2 ancestry, and loci SudRes1 and SudRes2 are additive effects per each copy of the minor allele at the representative codominant markers c582:65,596 and sc94:2,174,117. ***p < 0.001.

Fig. 3.

Characterization of Biomphalaria sudanica SudRes1 and SudRes2 genomic regions. (A) Pooled-GWAS p values in SudRes1 and SudRes2 regions (dashed boxes), which contain pooled-GWAS dual-variants (red) and other variants (blue), defined as in Figure 1. Contigs (grey rectangles) on B. sudanica chromosomes 5 and 6 are arranged horizontally based on contig orthology to 18 chromosomes of the B. glabrata genome (xgBioGlab47.1, NCBI RefSeq: GCF_947242115.1) and linkage map analysis (Dataset S2, Fig. S2–S4). Gene positions are shown (red/yellow/brown/orange boxes), highlighting three particularly prevalent classes of genes: in SudRes1 the multiple epidermal growth factor (MEGF) and galactose-binding like domain (GBD) containing receptor-like tyrosine-specific protein phosphatase (RPTP), other protein coding genes containing MEGF domains, and; in SudRes2 encoding a class of leucine-rich repeat containing G-protein couple receptors (GRL101) with C-type lectin and low-density lipoprotein extracellular domains (partial GRL101 genes included). (B) The predicted protein structure of receptor-like tyrosine-specific protein phosphatase (RPTP) coding gene BSUD.17727 (c6844) present in the SudRes1 region of the B. sudanica genome and containing intronic validated-variants (Fig. S8A). A similar RPTP coding gene is contained within adjacent contig c582 within SudRes1, and contig c1041 neighboring SudRes1 region (Fig. 3A). (C) A representative predicted protein structure of a GRL101-like G-protein coupled receptor coding gene, twelve of which were predicted through manual annotation within the SudRes2 region of contig sc94 (1.82–2.26 Mb) in the B. sudanica genome.