*2.5. Cone Bioassays*

Mosquitoes were exposed to PermaNet® 3.0 LN (Vestergaard Frandsen SA, Denmark) a LLIN consisting of a top panel made of monofilament polyethylene (100 denier) fabric incorporating deltamethrin at 4 g/kg (approx. 180 mg/m2) and piperonyl butoxide at 25 g/kg (approx. 1.1 g/m2), plus side panels made of multifilament polyester (75 denier) fabric with a strengthened border treated with deltamethrin at 2.8 g/kg (approx. 118 mg/m2) in WHO cone bioassays [25]. Following net airing of 2 weeks, pieces of netting (25 cm × 25 cm) were cut from the roof and side of the PermaNet 3.0 and cohorts of approximately 50 mosquitoes of each strain were exposed using the WHO standard protocol. Controls were exposed to insecticide free net in two replicates, each with 5 mosquitoes, one just before and one just after the treated exposures. Following exposure, the mosquitoes were aspirated into paper cups and supplied with 10% sucrose solution, and the initial knockdown effect was scored at 1 h and mortality was scored at 24 h post exposure.

#### *2.6. Target Site Mutation Genotyping*

Genomic DNA was collected within the first 5 months of colonisation of each strain and every subsequent 6–12 months thereafter. The DNA was extracted from 48 non-bloodfed females using a Qiagen blood and tissue DNA extraction kit (Qiagen, Germantown, MD, USA). Species ID was identified using the SINE PCR protocol [21].

Each strain was genotyped to identify the frequency of known target site resistance alleles (alleles 995F, 995S and 1570Y in the VGSC, the *ace-1*119S allele and the RDL alleles 296G and 296S) using Taqman™ assays [26–29]. The allelic variant 114T of the glutathione transferase *GSTe2* gene was also genotyped as previously [30].

#### *2.7. RNAseq Transcriptomic Analysis*

RNA was extracted from pools of 5, 3–5 day old presumed-mated adult females, snap frozen in the −80 ◦C at 10 am, using a PicoPure kit (Applied Biosystems Thermo Fisher, Waltham, MA, USA, after homogenisation with a motorised pestle. Quality and quantity of the RNA was analysed using an Agilent TapeStation (Santa Clara, CA, USA) and Nanodrop (Thermo Fisher) respectively, and three (Moz, Gaoua-ara, N'Gousso, Tiefora) or four (Banfora, VK72014, Kisumu, Bakaridjan) replicates from each strain sent for sequencing at Centre for Genomics, Liverpool, UK (RNA extractions for Banfora were performed as part of a separate study [31] but using the same methodology).

The resulting data was run through appropriate QC using FastQC and aligned to the latest *Anopheles gambiae* s.l. genome assembly PEST4 using Hisat2 with default parameters. The resulting bam file was sorted using samtools and the number of reads aligned to each gene extracted using featureCounts. Over 70% read assignment was seen for each replicate of each population with the majority showing >85%. Data from the *An. gambiae* s.s and *An. coluzzii* resistant populations were compared to the two susceptible populations (Kisumu and N'Gousso) using limma. First, a model matrix was defined to account for the populations and then contrasts were made to compare the resistant *An. gambiae* and *An. coluzzii* to both susceptible populations through the function makeContrasts using resistant—(N'Gousso + Kisumu)/2. Counts were then transformed to log2 counts per million reads (CPM), residuals calculated, and a smoothed curve fitted using the voom function which utilises normalisation factors calculated using calcNormFactors. lmFit was used to fit a linear model for each gene and eBayes used to smooth the standard errors. The function topTable was then used to retrieve results and written out to file; significance was taken as adjusted *p* value ≤ 0.05. In the case of the single *An. arabiensis* population, the contrast matrix was simply a resistant vs. susceptible design. In each instance the filter-ByExpr function from the EdgeR package was used to remove genes with low read number. Enrichment analysis was performed using the built-in GO term enrichment analysis on VectorBase with a Benjamini significance cut-off of ≤ 0.05. Revigo was then used to remove redundant GO terms allowing more appropriate visualisations; default parameters were used with a 0.5 selection. A custom table was also used with hypergeometric tests with fdr cut-off of ≤0.05 to integrate KEGG, Reactome and a priori genes of interest into the enrichment analysis (https://github.com/VictoriaIngham/BurkinsStrains) (accessed on 9 December 2021). All RNAseq data is deposited in SRA under accession PRJNA780362 and PRJNA750256.

#### *2.8. Metabolic Resistance—Detox Gene Expression Levels*

One to four μg of RNA, extracted from three pools of 5, 3–5-day-old female as described above, was reverse transcribed using Oligo dT (Invitrogen, Warrington, UK) and Superscript III (Invitrogen). The resulting cDNA was diluted to 4 ng/μ<sup>L</sup> and used as a template in the subsequent PCR reactions. Primers and probes as described by Maviridis et al. [32] were ordered from Integrated DNA Technologies (Leuven, Belgium), with Cy5 replacing Atto647N. Primers and probes were diluted to 10 μM for use in a 10 μL final reaction. Four multiplex reactions were carried out on each cDNA set in technical triplicate, as follows: (i) CYP6P4, CYP6Z1 and RPS7; (ii) CYP4G16 and CYP9K1; (iii) CYP6M2 and CYP6P1;

(iv) CYP6P3 and GSTE2. PrimeTime Gene Expression Master Mix (Integrated DNA Technologies) was used to set up each reaction following the manufacturer's instructions. Each reaction was carried out on a MxPro 3005P qPCR System (Agilent) with the following thermocycling conditions: 3 min at 95 ◦C followed by 40 cycles of 15 s at 95 ◦C; 1 min at 60 ◦C. Cycle threshold (Cq) values were exported and analysed using the ΔΔct methodology [33], using RPS7 as an endogenous control. Gaoua-ara were normalised against the susceptible Moz strain of *An. arabiensis*, and Bakaridjan, Banfora M, Tiefora and VK7 2014 were compared to both N'Gousso and Kisumu. A homogeneity of variance test was used to determine if data were normally distributed. Δct values were transformed to normalise (where applicable) and an ANOVA test, followed by Dunnett's test was performed. Where transformations did not normalise the data, a Dunn test was performed.
