Overexpression and nonsynonymous mutations of UDP-glycosyltransferases are potentially associated with pyrethroid resistance in Anopheles funestus

Abstract

UDP-glycosyltransferases (UGTs) enzymes are pivotal in insecticide resistance by transforming hydrophobic substrates into more hydrophilic forms for efficient cell elimination. This study provides the first comprehensive investigation of Anopheles funestus UGT genes, their evolution, and their association with pyrethroid resistance. We employed a genome-wide association study using pooled sequencing (GWAS-PoolSeq) and transcriptomics on pyrethroid-resistant An. funestus, along with deep-targeted sequencing of UGTs in 80 mosquitoes Africa-wide. UGT310B2 was consistently overexpressed Africa-wide and significant gene-wise F_st differentiation was observed between resistant and susceptible populations: UGT301C2 and UGT302A3 in Malawi, and UGT306C2 in Uganda. Additionally, nonsynonymous mutations in UGT genes were identified. Gene-wise Tajima's D density curves provide insights into population structures within populations across these countries, supporting previous observations. These findings have important implications for current An. funestus control strategies facilitating the prediction of cross-resistance to other UGT-metabolised polar insecticides, thereby guiding more effective and targeted insecticide resistance management efforts.

Keywords: Genomics; Insecticide resistance; Target sequencing; Transcriptomics; UDP-glycosyltransferases; Vector control.

Fig. 1

Phylogenetic tree of 123 UGT sequences from An. funestus, Ae. Aegypti, An. gambiae, D. melanogaster. Phylogenetic analysis was constructed using FastTree method. Fast tree support values that are larger than 0.5 are marked on the corresponding branch. UGT names start with the species initials An. funestus (Af), Ae. Aegypti (Aa), An. gambiae (Ag), and D. melanogaster (Dm), followed by the family classification (numbers), the subfamily classification (letters), and the unique numbers. The numbers on the outer circle refer to the family classification.

Fig. 2

Differentially expressed An. funestus UGT genes in response to permethrin resistance. (A) Differentially expressed genes (DEGs) between the field-resistant population after exposure to permethrin and the unexposed population from the same country were identified in Malawi but not in Cameroon and Uganda. (A) Comparing the transcriptional profile of resistant populations from Malawi, Cameroon, Uganda and FUMOZ to that of susceptible FANG population identified Africa-wide upregulation of UGT310B2. The red arrow pointing upwards indicates upregulation and the blue arrow pointing downwards indicates downregulation. (B) All An. funestus UGTs (x-axes) plotted against fold change (FC) (y-axes) obtained from DEGs analyses colour coded according to the plot legend highlights the Africa-wide upregulation of UGT310B2 with a maximum fold change of 5.6 in resistant population from Malawi when compared with transcription profile of FANG. UGT310B2, UGT306D2, UGT301E3, and UGT315A3 have a significant overexpression in Malawi resistant population when compared with the FANG profile and with the unexposed population from Malawi highlighting the potential role of UGTs in detoxification in Malawian population. (C) The overexpression of UGT310B2 was detected in the FUMOZ colony using qPCR. However, due to the relative overexpression of UGT301C2 from the FUMOZ RNAseq profile compared to the FANG detecting the overexpression in the FUMOZ colony by qPCR was challenging. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 3

Gene-wise F_ST for all genes included in the targeted sequencing in all analyses between resistant and susceptible. The average and 0.95 quantiles of gene-wise F_ST for each chromosome are represented by the blue horizontal dotted line and the red dotted line respectively. The genomic location of all genes in the targeted region (x-axis), represented by a circle, was plotted against the Gene-wise F_ST (y-axis). UGT genes included in the targeted sequencing are colour coded according to the plot legend. Raw data of all the genes included in the targeted sequencing is in (Supplementary Dataset 5), averages and 0.9 quantiles for each comparison are in (Table S2), and UGT genes F_ST values are in (Table S3). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 4

Tajima's D for UGT genes and density curves of all the genes included by the targeted enrichment sequencing for all the populations. In populations from Southern Africa Malawi (MAL) (orange), FANG (Angola) (Yellow) and FUMOZ (Mozambique) (Blue) Tajima's D density curves were close to equilibrium (Black) (B). In Malawi, UGT301C2 is impacted by a recent selective sweep with a negative below the 0.05 quantiles of simulated Tajima's D (A). In Uganda (UG) and Cameroon (CAM), gene-wise Tajima's D was predominantly skewed towards negative values. In Cameroon, 7 UGT genes are below 0.05 quantiles of simulated Tajima's D and most UGT genes in Uganda are with a negative Tajima's D but not within the lowest 5% (A). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Fig. 5

Targeted UGTs polymorphism and haplotype network analysis. UGTs gene-wise nucleotide diversity was calculated between (resistant and susceptible) populations and within populations from Malawi, Cameroon, Uganda and laboratory colonies the FANG and the FUMOZ (A). The number of polymorphic substitutions for UGTs in each country reflects the gene size (Table 1) (B). TCS haplotype network reveals limited haplotype clustering in UGT301C2 (C), UGT306A3 (D) and UGT306C2 (E) where predominant haplotypes shared by different populations were detected. The predominant haplotype in UGT301C2 was observed in (40/140) haplotypes from Cameroon (9), Uganda (14), Malawi (12) and FUMOZ (5). While 18/20 FANG haplotypes cluster together in three predominant nodes separated by a single mutation step each contains (9/20), (7/20), and (2/20) haplotypes respectively (Fig. 5C). In UGT306A3 two main haplotype clusters separated by a single mutation step were detected. The biggest cluster contains (45/160) haplotypes where all the 20 FANG haplotypes cluster, 9/20 FUMOZ and 15/40 haplotypes from Malawi. The second cluster contains (41/160) haplotypes including 12/40 from Malawi, 16/40 from Cameroon and 7/40 from Ugandan population (Fig. 5D). In UGT306C2 many haplotype clusters were observed the most dominant is 23 haplotypes cluster specific to FUMOZ (7/20), FANG (6/20) and Malawi (11/40).

Fig. 6

Significantly differentiated UGT genes nonsynonymous SNPs between laboratory colonies and in each country. Non-synonymous SNPs that are significantly differentiated between susceptible and resistant in each analysis are colour-coded according to the associated UGT gene. The red dotted line indicates the significance level at p-value 0.05 (−log(0.05) = 2.99). Further details of those nonsynonymous SNPs are included in Table S8 and the location of those SNPs relative to UGTs conserved sites is highlighted in Supplementary Fig. 11. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)