The Pattern of Knockdown Resistance Mutations Correlated to the Annual Average Temperature in Field Populations of Aedes Albopictus (Diptera: Culicidae) in China

Background: Aedes albopictus is the main vector of dengue fever in China, distributed from north to south in China. Insecticides are an important method to control the mosquitoes, especially in the outbreak of dengue fever, but insecticide resistance raises the risk of failure to control vector-borne diseases. Knockdown resistance (kdr) caused by point mutations in the VGSC gene is a key mechanism that confers resistance to pyrethroids. To explore the characteristics and possible evolution trend of kdr mutation in Ae. albopictus, we analyzed the kdr mutations of eld populations in China in this study. Methods: A total of 1 549 Ae. albopictus were collected from 18 sites in China from 2017 to 2019, as well as 50 individuals from three sites in the 1990s. A fragment of approximately 350 bp from part of S6 segment in the VGSC gene domain III was amplied and sequenced. The haplotypes of VGSC gene were recorded and the parsimony network was constructed using TCS 1.21. The data of annual average temperatures (AAT) of collection sites was acquired from national database. The correlation between AAT of the collection site and the kdr mutation rate was analyzed by Pearson Correlation using SPSS 21.0. The overall of 45.62%. Nine mutant alleles were detected at codon 1534 in fteen eld populations, namely TCC/TCG and Only one in multiple mutations at I1532 and in a sample mutation signicantly positive AAT (Coecient=0.624, p=0.0056), the 1532 mutation rate signicantly to AAT (Coecient=-0.645, p=0.0038). mutation


Background
Dengue fever is one of the most widespread mosquito-borne diseases around the world. It was estimated that there are 390 million dengue virus infections worldwide each year (1). Dengue fever has also been an important public health emergency in China. In 2019, 22 599 dengue fever cases were reported in China, and the distribution area reached an unprecedented level (2,3). Due to the scarcity of therapeutic drug and effective vaccine, the vector controlling is still the main measure to prevent dengue fever (4).
Aedes aegypti and Ae. albopictus are the main vectors of dengue fever (5). Although Ae. aegypti is more susceptible to dengue virus, Ae. albopictus is considered the main vector of dengue fever in China due to its wide distribution range (6). Ae. albopictus is distributed from Liaoning Province in the north to Hainan Province in the south in China, covering tropical, subtropical and temperate climate (6). Insecticide is an essential method to control the Ae. albopictus vector, especially during outbreaks of mosquito-borne diseases.
With the advantages of high e ciency and low toxicity, pyrethroid insecticides have been widely used to control mosquito in China since 1980s (7). Although integrated control are emphasized in the general Ae. albopictus control, pyrethroid insecticides were sprayed heavily during the outbreak of dengue fever to control the epidemic quickly(8, 9). Long-term and mass application of insecticides has led to a serious problem --insecticide resistance. High level and multiple pyrethroid insecticides resistance have been observed in many Ae. albopictus eld populations in China (10)(11)(12).
The target of pyrethroids is voltage-gated sodium channel (VGSC), which interferes with electronic signals in the nervous system, leading to paralysis and death, an effect known as knockdown (13). The target insensitivity was a critical mechanism of the pyrethroids resistance developing, which was caused by mutations in the VGSC, so-called knockdown resistance mutation (kdr mutation) (14,15).
More and more kdr mutations have been discovered, some of which have been con rmed to be closely related to insecticide resistance (12,31,32). However, how these mutations appear, accumulate and spread in eld populations remain a puzzling and attractive issue. Moreover, the evolution trend of kdr mutation maybe the key which would determine the strategy of pesticide application in the future.
In this study, we collected eighteen eld populations from North to South China, as well as three populations in the 1990s. The types and frequencies of kdr mutations were surveyed and analysis in order to explore the characteristics and possible evolution trend of kdr mutation.

Sample collection and species identi cation
Eighteen eld populations of Ae. albopictus were collected from in 13 provinces (municipalities) in China from 2017 to 2019 (Fig. 1). Larvae and pupae were scooped from breeding sites and brought back to laboratory, which were reared to adults under standard conditions at 26 ± 1°C and 65 ± 5% relative humidity with a photoperiod of 12-h light: 12-h dark. Adult mosquitoes were collected using sucking aspirators, BG-trap (Biogents, Germany) or light trap (Houji Dianzi, China) in outdoor environment near human. To explore the historical changes of kdr mutations in Ae. albopictus, three eld populations in the 1990s were also collected, including HNHK94 (Haikou, 1994), GDGZ95 (Guangzhou, 1995) and SCCD97 (Chengdu, 1997). The collection information was summarized in Table 1. Species of Aedes adult mosquitoes were identi ed by morphology characteristics (33) and con rmed by molecular markers (34). *These population's information has been published in our previous work [12,30]. & The eld populations collected in 1990s.

Annual average temperature of collection sites
The annual average temperature (AAT) of collection sites were search and recorded from China National Meteorological Science Data Center (http://data.cma.cn) (Table S1).

DNA extraction
Genomic DNA was extracted from adult mosquitoes using DNAzol Reagent (Invitrogen, USA) according to the manufacturer's instructions. In brief, a single mosquito was put in a 1.5mL eppendorf tube, and powdered e ciently by rstly freezing in liquid nitrogen. Then 100µL DNAzol were added to lysis, and the supernatant was moved to another eppendorf tube. After adding 50µL ethanol, the mixed solution was stored at -20℃ overnight. The DNA in the solution was adsorbed through a nucleic acid adsorption column. After washing with 70% ethanol, the DNA was eluted from the column with ddH 2 O. The concentration of the DNA was quanti ed using OD 280 absorption.

PCR and kdr alleles detection
A fragment of approximately 350 bp of S6 segment in the VGSC gene domain III was ampli ed using PCR kit (Aidlab, China) and sequenced. The primers are aegSCF7 (5'-AGG TAT CCG AAC GTT GCT GT-3') and aegSCR8 (5'-TAG CTT TCA GCG GCT TCT TC-3') (12,16). PCR reaction was carried out in a Verity 96 well 157 Thermal Cycler (Applied Biosystems, USA). The cycling parameter included an initial step of denaturation at 94°C for 2 min, followed by 35 cycles of ampli cation at 94℃ for 30 s, 52℃ for 30 s, and 72℃ for 30 s, followed by a nal extension step at 72℃ for 8 min. After electrophoresis, PCR products were puri ed and sequenced from both directions. Sequences were aligned and analyzed by DNASTAR Lasergene 12.0 software (35).

Parsimony network the haplotypes of VGSC gene
The haplotypes of VGSC gene were recorded and the networks was constructed based on statistical parsimony using TCS 1.21 software (36).

Statistical analyses
The allele frequency of kdr mutations was calculated as: "R" means the allele frequency. "k" represents the number of alleles, and "n" represents the sample size.
The mutation frequency was de ned as the frequency of wildtype-mutant heterozygotes and mutant-mutant homozygotes, which was calculated as: "M" means the mutation frequency. "a" represents the number of wildtype-mutant samples, and "b" represents the number of mutant-mutant samples. "n" is the sample size.
The correlation between the AAT of the collection site and the kdr mutation rate was analyzed by Pearson Correlation using SPSS 21.0. The difference of the mutation rate between groups was calculated by one-way ANOVA using SPSS 21.0.
In terms of the mutant allele types, the SZSK population carried the most types with eight mutant alleles, followed by GXNN, YNJH and HNSY with 4 mutant alleles, GDGZ and SHYP with 3 mutant alleles, HNHK and SDJNan with 2 mutant alleles. For the other seven populations, only one mutant allele as F1534S/TCC was detected ( Table 2).
1534TCC/S is the most common mutant allele, which was detected in all of the 15 populations. A total of four mutant alleles were found in only one population, namely CTG/L, CGC/R, TGG/W and TTA/L. Among them, CTG/L (3.33%), CGC/R (6.67%) and TGG/W (2.50%) were only found in SZSK, and TTA/L (10.65%) was identi ed in the SDJNan.

Multiple mutations at codon 1532 and 1534
Multiple mutations at I1532 and F1534 appeared in ve populations as YNJH, SHBS, SHYP, SHGQ and ZJHZ (Fig S3). The highest frequency of cooccurrence of the mutations was 21.43% in SHYP population, followed by 15.08% (SHGQ). Three types of combined mutations were present, as I/T + F/S, I/T + S/S, T/T + F/L. However, no individuals with mutant homozygotes at both codons were found (Table 3).

Geographic distribution of kdr mutants of Ae. albopictus in China
Interestingly, the pattern of kdr mutations in Ae. albopictus seems to be related to geographical factors. F1534S/L/C/R/W mutant alleles were mainly found from southern and central China ( Fig. 2A), while I1532T was mainly found in populations from northern and central China (Fig. 2B). If the eld populations were divided to the southern group and the northern group by 32°N (Fig. 2), there is a signi cant difference in the 1534 mutation rate (p = 0.0003) and the 1532 mutation rate (p = 0.0443) between the two groups.
To further explore the phenomenon, we collected the AAT of the collection sites. The 1534 mutation rate was signi cantly positive related to AAT (Coe cient = 0.624, p = 0.0056), while the 1532 mutation rate was signi cantly negative related to AAT (Coe cient=-0.645, p = 0.0038).
It is worth noting that the ZJHZ population was collected during the 2019 Hangzhou dengue fever epidemic. At that time, a large number of insecticides had been used, so the collected mosquitoes were likely to be a screened resistant population. The collection site of the GDST population is very close to the sea. It was learned that almost no pesticides was used in the local area through oral interviews. If the ZJHZ and GDST population were removed from the correlation analysis, the 1534 mutation rate showed more signi cantly positive related to AAT (Coe cient = 0.843, p < 0.0001), and the 1532 mutation rate was also signi cantly negative related to AAT (Coe cient=-0.657, p = 0.0057).

Genealogical relationships among the haplotypes of VGSC gene
Based on the alleles in codons 1532 and 1534, thirteen haplotypes were inferred, which has been submitted to the NCBI (https://www.ncbi.nlm.nih.gov/) ( Table 4). The 95% parsimony network showed a parallel evolution pathway. The ancestral wildtype haplotype (H01_IF) weight was only 0.13. Six mutant haplotypes were formed by one step mutation from ancestral haplotype, as H02_IL, H04_IL, H05_IL, H06_IS, H08_IC and H11_TF. The other six haplotypes were formed by one more mutation (Fig. 3). This result supports the view that kdr mutations of Ae. albopictus have multiple origins.

Historic change of kdr mutants in Ae. albopictus in China
A total of 50 individuals from three sites in the 1990s were collected to explore the historic change of the kdr mutation, including HNHK94 (Haikou, 1994, n = 12), GDGZ95 (Guangzhou, 1995, n = 18) and SCCD97 (Chengdu, 1997, n = 20)). No mutant allele was detected at codon 1532 of VGSC gene. However, F1534S/TCC was found in HNHK94 with an unexpected frequency of 100%. HNHK94 samples were collected from Haikou in 1994, when15 years earlier than the rstly reported kdr mutant populations by Kasai et al (16). Although the number of samples available is limited, the result clearly showed that the F1534S kdr mutation of Ae. albopictus appeared no later than 1994.

Synonymous mutations
In addition to above non-synonymous mutations, we also found 5  (Table S4). However, we did not nd any correlation between these synonymous mutations and non-synonymous mutations at codon 1534 and 1532. The role of these synonymous mutations remains to be explored.

Discussion
Insecticide resistance has become an important threat to the control of diseases transmitted by Ae. albopictus. In this study, we collected the eld populations of Ae. albopictus in China, detected and analyzed the kdr mutation to explore the characteristics and possible evolution trend of the kdr mutations. The sample collection sites covered the main distribution range of Ae. albopictus in China from 16°N to 40°N in this study (6). The climate includes tropical, subtropical and temperate climates.
The results showed that types of kdr mutations have appeared widely in the eld populations of Ae. albopictus in China. Certain population (SZSK) showed a high mutation rate of 1534 with complex mutation types. Moreover, two novel mutant alleles F1534W/TGG and F1534R/CGC were detected in SZSK population. The AAT of this group collection point is relatively high (23.0°C), where mosquitoes reproduce all year round. Moreover, developed trade and logistics, high population mobility, and frequently used pesticides in Shenzhen are all factors that promote insecticide resistance in mosquito populations (37). If the use of pesticides cannot be controlled, the high frequency and multiple types of 1534 mutations appeared in SZSK population may be the trend of kdr mutations in the eld populations of Ae. albopictus in the future.
An extremely high mutation rate of 1534 was detected in the HNHK94 population collected in Haikou City in 1994, indicating that some populations of Ae. albopictus had high insecticide resistance decades ago. The result pushed the discovery time of kdr mutation in Ae. albopictus samples forward by 15 years (16). The phenomenon may be related to the large-scale use of pesticides in Hainan Province after three dengue fever outbreaks in the last century, as 1979 1982, 1985 ~ 1988, 1991(38). The ZJHZ population was also collected after the occurrence of dengue fever in Hangzhou (11), which similarly showed a very high percentage of 1534 mutation rates, suggesting that the large-scale use of pesticides is a "closely related" factor for the 1534 mutation rate of the Ae. albopictus population. Therefore, a more cautious, precise and strict insecticide resistance management (IRM) is urgently to be established (39).
In northern populations where the average temperature was low, the kdr mutation rate of Ae. albopictus was generally low, and the I1532 mutation rate was relatively high. However, in the southern populations with higher average temperature, the kdr mutation rate was higher, and the F1534 mutation rate was relatively high. Research on the correlation between gene mutations and insecticide resistance showed that F1534S and F1534C were signi cantly related to insecticide resistance (12,21), while the correlation between I1523T and insecticides was still unclear (12,32). The results of this study clearly show that the southern population is more resistant to pesticides than the northern population. The temperature in the south is high, and mosquitoes have more generations, which can be regarded as fast-evolving populations, while the northern populations can be regarded as slow-evolving populations. The results in this study suggested that the highly resistant population may be screened out in a relatively short period of time in the fast-evolving populations. Due to the long cold winter in the north, pesticides are often used seasonally, which may slow down the development of pesticide resistance. This suggests that seasonal use of pesticides may be a strategy to slow down the development of pesticide resistance. However, how to use insecticides seasonally in a certain area and effectively control the risk of vector-borne disease transmission requires further research. This is very important for the prevention and control of mosquito-borne diseases in the context of global climate change (40).
Although some meaningful results were obtained, there are still some limitations in this study. First, there are only a few historical samples available, which make it impossible to systematically analyze the historical changes of kdr mutations. In addition, accurate information on the use of insecticides cannot be collected, which affected the interpretation of the results. Moreover, this study focused on the most common kdr mutations in Ae. albopictus, leaving some types of kdr mutations with low mutation rate (18). Finally, although the synonymous mutations and intron length polymorphisms have been documented in the populations of Ae. albopictus and Ae. aegypti (19,28,41), the signi cance of the synonymous mutations found in this study is not clear, and further research is needed.

Conclusions
From this study, we can found the kdr mutations are widespread in the eld populations of Ae. albopictus in China. Two novel mutant alleles F1534W/TGG and F1534R/CGC were the detected for the rst time. The 1534 kdr mutation appeared in the population of Ae. albopictus no later than 1994. HNHK94 was the earliest known Ae. albopictus eld population with kdr mutation. F1534 mutation rate is positive correlated to AAT, while I1532 mutation rate is negative correlated to AAT. Insecticide using should be carefully managed to slow down the spread of high-resistance Ae. albopictus populations. Availability of data and material

Abbreviations
The datasets used and/or analyzed during the current study are available from the corresponding authors on reasonable request, and sequences are available in GenBank.

Figure 1
A schematic map of collection sites for Aedes albopictus in China and the composition of the codon 1534 alleles of VGSC gene in the samples. The red balloons indicate the sampling sites (see Table 1 for details). The codon 1534 alleles composition in each site was displayed by a pie chart with green as wildtype F1534, yellow as mutant allele F1534S and blue as mutant allele F1534L/C/R/W (see Table 2 for details). Note: The designations employed and the presentation of the material on this map do not imply the expression of any opinion whatsoever on the part of Research Square concerning the legal status of any country, territory, city or area or of its authorities, or concerning the delimitation of its frontiers or boundaries. This map has been provided by the authors.

Figure 2
The histogram of allele frequency at codon 1532 (A) and 1534 (B) of VGSC gene in Aedes albopictus populations. SH was the pool of three populations (SHBS, SHYP, SHGQ). For information details, see Table 1 and Table 2. The populations on the left of the dotted line were from central and northern China, and the populations on the right side were from southern China. Horizontal axis indicates the populations while vertical axis denotes the allele frequency. The populations with asterisk was co-occurrence of mutations both codons 1532 and 1534.