A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth

A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth | bioRxiv /* */ /* */ <!-- <!-- /*! * yepnope1.5.4 * (c) WTFPL, GPLv2 */ (function(a,b,c){function d(a){return"[object Function]"==o.call(a)}function e(a){return"string"==typeof a}function f(){}function g(a){return!a||"loaded"==a||"complete"==a||"uninitialized"==a}function h(){var a=p.shift();q=1,a?a.t?m(function(){("c"==a.t?B.injectCss:B.injectJs)(a.s,0,a.a,a.x,a.e,1)},0):(a(),h()):q=0}function i(a,c,d,e,f,i,j){function k(b){if(!o&&g(l.readyState)&&(u.r=o=1,!q&&h(),l.onload=l.onreadystatechange=null,b)){"img"!=a&&m(function(){t.removeChild(l)},50);for(var d in y[c])y[c].hasOwnProperty(d)&&y[c][d].onload()}}var j=j||B.errorTimeout,l=b.createElement(a),o=0,r=0,u={t:d,s:c,e:f,a:i,x:j};1===y[c]&&(r=1,y[c]=[]),"object"==a?l.data=c:(l.src=c,l.type=a),l.width=l.height="0",l.onerror=l.onload=l.onreadystatechange=function(){k.call(this,r)},p.splice(e,0,u),"img"!=a&&(r||2===y[c]?(t.insertBefore(l,s?null:n),m(k,j)):y[c].push(l))}function j(a,b,c,d,f){return q=0,b=b||"j",e(a)?i("c"==b?v:u,a,b,this.i++,c,d,f):(p.splice(this.i++,0,a),1==p.length&&h()),this}function k(){var a=B;return a.loader={load:j,i:0},a}var l=b.documentElement,m=a.setTimeout,n=b.getElementsByTagName("script")[0],o={}.toString,p=[],q=0,r="MozAppearance"in l.style,s=r&&!!b.createRange().compareNode,t=s?l:n.parentNode,l=a.opera&&"[object Opera]"==o.call(a.opera),l=!!b.attachEvent&&!l,u=r?"object":l?"script":"img",v=l?"script":u,w=Array.isArray||function(a){return"[object Array]"==o.call(a)},x=[],y={},z={timeout:function(a,b){return b.length&&(a.timeout=b[0]),a}},A,B;B=function(a){function b(a){var a=a.split("!"),b=x.length,c=a.pop(),d=a.length,c={url:c,origUrl:c,prefixes:a},e,f,g;for(f=0;f<d;f++)g=a[f].split("="),(e=z[g.shift()])&&(c=e(c,g));for(f=0;f<b;f++)c=x[f](c);return c}function g(a,e,f,g,h){var i=b(a),j=i.autoCallback;i.url.split(".").pop().split("?").shift(),i.bypass||(e&&(e=d(e)?e:e[a]||e[g]||e[a.split("/").pop().split("?")[0]]),i.instead?i.instead(a,e,f,g,h):(y[i.url]?i.noexec=!0:y[i.url]=1,f.load(i.url,i.forceCSS||!i.forceJS&&"css"==i.url.split(".").pop().split("?").shift()?"c":c,i.noexec,i.attrs,i.timeout),(d(e)||d(j))&&f.load(function(){k(),e&&e(i.origUrl,h,g),j&&j(i.origUrl,h,g),y[i.url]=2})))}function h(a,b){function c(a,c){if(a){if(e(a))c||(j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}),g(a,j,b,0,h);else if(Object(a)===a)for(n in m=function(){var b=0,c;for(c in a)a.hasOwnProperty(c)&&b++;return b}(),a)a.hasOwnProperty(n)&&(!c&&!--m&&(d(j)?j=function(){var a=[].slice.call(arguments);k.apply(this,a),l()}:j[n]=function(a){return function(){var b=[].slice.call(arguments);a&&a.apply(this,b),l()}}(k[n])),g(a[n],j,b,n,h))}else!c&&l()}var h=!!a.test,i=a.load||a.both,j=a.callback||f,k=j,l=a.complete||f,m,n;c(h?a.yep:a.nope,!!i),i&&c(i)}var i,j,l=this.yepnope.loader;if(e(a))g(a,0,l,0);else if(w(a))for(i=0;i (function(w,d,s,l,i){w[l]=w[l]||[];w[l].push({'gtm.start':new Date().getTime(),event:'gtm.js'});var f=d.getElementsByTagName(s)[0];var j=d.createElement(s);var dl=l!='dataLayer'?'&l='+l:'';j.src='//www.googletagmanager.com/gtm.js?id='+i+dl;j.type='text/javascript';j.async=true;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-M677548'); Skip to main content Home About Submit ALERTS / RSS Search for this keyword Advanced Search New Results A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth View ORCID Profile Pierre Nouhaud , View ORCID Profile Mathieu Gautier , View ORCID Profile Alexandre Schifano , View ORCID Profile Laure Sauné , View ORCID Profile Carole Kerdelhué , View ORCID Profile Charles perrier doi: https://doi.org/10.1101/2025.11.24.690095 Pierre Nouhaud 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Pierre Nouhaud For correspondence: pierre.nouhaud{at}inrae.fr Mathieu Gautier 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Mathieu Gautier Alexandre Schifano 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France 2 Current address: Program in Evolutionary Biology, Department of Ecology and Genetics (IEG), Evolutionsbiologiskt centrum , Norbyvägen 18D, Uppsala University, 752 36 Uppsala, Sweden Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Alexandre Schifano Laure Sauné 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Laure Sauné Carole Kerdelhué 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Carole Kerdelhué Charles perrier 1 UMR CBGP, INRAE, CIRAD, IRD, Institut Agro, Université Montpellier , Montpellier, France Find this author on Google Scholar Find this author on PubMed Search for this author on this site ORCID record for Charles perrier Abstract Full Text Info/History Metrics Supplementary material Preview PDF Abstract Understanding how natural populations adapt to complex environmental gradients is crucial for predicting evolutionary responses to global change. The pine processionary moth ( Thaumetopoea pityocampa ), a major forest pest expanding northward and upward in Europe, provides an ideal model to explore genomic adaptation along altitudinal and/or latitudinal gradients. We combined pooled and individual whole-genome resequencing of four pairs of low- and high-elevation populations from Spain, Italy and France (mainland and Corsica) to detect signatures of local adaptation. Population structure analyses revealed strong differentiation among distant populations but limited divergence within altitude pairs. Genome scans identified a few candidate regions under selection, and we notably uncovered a large (∼1 Mbp) region showing reduced nucleotide diversity, negative Tajima’s D , and fixed allele differences between high- and low-altitude populations, consistent with a recent selective sweep. This region includes a cluster of cytochrome P450 genes and the voltage-gated sodium channel gene para , both involved in detoxification and insecticide resistance. Although signatures of selection were observed at this locus in two distant population pairs, no definitive evidence of recent introgression was found between these distant populations, suggesting independent evolution. These results indicate that while climate variation may drive adaptation in T. pityocampa , local insecticide exposure may exert a stronger selective pressure and mask other genomic signals. Our study highlights candidate genes and regions for further investigation of resistance evolution in an insect species in the wild. Significance statement Understanding how species adapt to their local environments is critical for predicting their survival under changing conditions. Here we show that populations of the pine processionary moth exhibit a strong genomic signature of local adaptation to high- and low-elevation habitats at a region containing genes involved in detoxification and insecticide resistance. These findings suggest that local environmental pressures, such as insecticide exposure, can drive rapid evolutionary changes, highlighting the importance of considering local adaptation when managing pest species and predicting their responses to environmental change. Introduction Climate change exerts a strong pressure on the vast majority of life forms, its swift pace increasing decline rate and extinction risk in natural populations ( Thomas et al. 2004 ; Bellard et al. 2012 ; Scheffers et al. 2016 ). To persist, species must either shift their distributions latitudinally and/or altitudinally, acclimatize to novel stressful conditions through phenotypic plasticity, and/or undergo rapid adaptation ( Parmesan & Yohe 2003 ; Chen et al. 2011 ; Inouye 2022 ). These strategies are not mutually exclusive, and natural populations may often combine different responses to survive, for instance through local adaptation at expanding range margins ( Hill et al. 2011 ). Yet, climate is not the sole driver of change, as colonization of new habitats also exposes populations to additional selective pressures, many of which may also be anthropogenically induced (for a review in insects, see McCulloch & Waters 2023 ). Population genetics provides a powerful framework to map signatures of selection along the genome, offering important insights into how natural populations may respond to climate change ( Bernatchez et al. 2024 ). Numerous methods have been developed to identify candidate adaptive loci from sequencing data (reviewed in, e.g., Bourgeois & Warren 2021 ). Divergence-based genome scans, for instance, have been widely used to detect loci whose genetic differentiation exceeds expectations under neutrality (e.g., Dudaniec et al. 2018 ). Alternatively, genotype–environment association (GEA) analyses assess statistical relationships between allele frequencies and environmental variables to identify the drivers of local adaptation ( Forester et al. 2018 ). Because both divergence and GEA scans are prone to type-I errors (see, e.g., Lotterhos & Whitlock 2015 ), taking into account the demographic history of populations ( Bonhomme et al. 2010 ; Günther & Coop 2013 ) and incorporating replicated population pairs (e.g., Pélissié et al. 2022 ) altogether provide a robust strategy to more confidently detect candidate adaptive loci showing signatures of selection. The pine processionary moth Thaumetopoea pityocampa (Dennis & Schiffermüller; Lepidoptera: Notodontidae) is one of the most severe pests of pine trees in Europe and North Africa (Kanat et al. 2005). Late instar larvae also represent a significant public health concern due to their urticating hairs, which can cause serious dermatological, ocular, and respiratory lesions in humans and animals ( Vasseur et al. 2022 ; Rivière 2011 ). The potential pine processionary moth distribution is mainly constrained by winter temperature (but see Rossi et al. 2025 ) as larvae are unable to feed and eventually die of starvation below ca. -15°C ( Robinet et al. 2015 ). While the species has been expanding across Europe since the last glacial maximum ( Kerdelhué et al. 2009 ), in recent decades, rising temperatures have increased the expansion pace of populations northward across Europe ( Battisti et al. 2005 ), along with a retraction of its southern range (Tunisia, Bourougaaoui et al. 2021 ). Whereas the northern distribution limit remained near 48°N until the 1980s, it now extends beyond 49°N ( Robinet et al. 2025 ), and species distribution models suggest that the species could expand to 50°N under future climate scenarios ( Rossi et al. 2025 ). In addition to its northward expansion, the pine processionary moth has been recorded establishing at higher elevations, with altitudinal colonization rates estimated at approximately 70 meters per decade on south-facing slopes (29 meters per decade on north-facing slopes, Battisti et al. 2005 ). Although the latitudinal and altitudinal range expansion of the pine processionary moth across Europe is well documented, with, e.g., earlier emergence recorded at higher elevation levels ( Martin et al. 2022 ), whether populations locally adapt to ecological variations on altitudinal gradients remains an open question. Here, we take advantage of the recent release of a chromosome-level assembly ( Gautier et al. 2025 ) to study potential genomic footprints of local adaptation to altitudinal ecological variation. We implemented a replicated sampling design, collecting individuals from low- and high-elevation sites colonized post-glacial warming (i.e., prior to the last few decades) across four distinct localities around the Mediterranean basin ( Figure 1 , Table 1 ). To obtain genome-wide data for each of these eight populations, we used a combination of pooled sequencing (pool-seq, Futschik & Schlötterer 2010 ) and individual sequencing (ind-seq). We first characterized the genetic structure of pine processionary moth populations. We then applied model-based approaches to identify candidate loci associated with local adaptation and with altitude, hypothesizing that genes involved in the circadian rhythm and in metabolism regulation could represent promising candidates. While we did not detect convergent signals of association to altitude across the four altitudinal gradients, we uncovered evidence of a strong putative selective sweep in two population pairs, with gene content suggesting a role in insecticide resistance. Altogether, our results indicate that multiple selective pressures are acting in natural populations and highlight the need for careful interpretation of adaptive signals in genomic data. Download figure Open in new tab Figure 1. Sampling locations and genetic structure of the pine processionary moth individuals and pools. a. Map of the sampled localities. b. PCA of the 89 ind-seq libraries genotyped at 84,443 linkage-pruned autosomal SNPs. c. Genome-wide nucleotide diversity (π, ×10 3 , diagonal) and pairwise divergence metrics ( F ST , ×10 2 , upper triangle; d XY , ×10 3 , lower triangle) computed on ind-seq libraries. d. and e. Individual supervised clustering results obtained with Admixture for K = 4 and K = 7 ancestry components, respectively. f. Hierarchical clustering tree of the 8 pool-seq libraries based on the scaled covariance matrix of population allelic frequencies (Ω) estimated by BayPass. View this table: View inline View popup Download powerpoint Table 1. Overview of sampling locations and associated data sets. Pool sizes are indicated as diploid, and numbers in brackets for ind-seq samples correspond to sex counts inferred from the Z-chromosome-to-autosomal depth ratio analysis (see Methods). Results Strong structure among distant populations but little differentiation within low and high altitude pairs Sequencing produced between 92.3-682 (median: 123) and 260-577 (median: 374) million reads for ind-seq and pool-seq libraries, respectively, with mean insert sizes ranging from 340 to 457 (median: 387, see Suppl. Table 1 for detailed sequencing statistics). Mean depths ranged from 15.7× to 109× (median: 21.9×) for ind-seq libraries and from 41.4× to 97.7× (median: 60.4×) for pool-seq libraries. Heterosomal-to-autosomal depth ratios identified 66 males and 23 females from the ind-seq data ( Table 1 , Suppl. Table 1, Suppl. Fig. 1). Joint SNP calling and subsequent library-specific filtering recovered 14,323,239 sites for the pool-seq dataset and a subset of 9,193,257 sites for the ind-seq dataset. The PCA carried on ind-seq data using a subset of 84,443 linkage-pruned autosomal SNPs clustered samples by locality, with the first axis discriminating between Corsican samples and the rest ( Fig. 1b ). Nucleotide diversity measured with π was highest in the two Spanish populations, intermediate in Mont Ventoux and comparably smaller in Corsica and Italy ( Fig. 1c ). Overall, absolute and relative pairwise divergence metrics were low to moderate between populations within localities, and the increased differentiation observed for the Italian pair could be caused by the greater distance between sampling sites (ca. 135 km, Fig. 1a ). Divergence was however higher between localities, with F ST values ranging from 0.11 (Mont Ventoux vs . Spain) to 0.38 (Corsica vs . Italy). Supervised clustering grouped individuals per locality at K = 4 ( Fig. 1d ), while Puechmaille’s method gave K = 7 as the most likely number of ancestry components ( Fig. 1e ). Under this scenario, both Corsican and Italian samples were split according to elevation levels, while Mont Ventoux samples remained within a single cluster and Spanish samples appeared as admixed, which could be a consequence of the strong nucleotide diversity present in these samples compared to the rest ( Fig. 1c ). The scaled covariance matrix of population allelic frequencies Ω estimated from pool-seq data using BayPass recapitulated the population structure identified from ind-seq data, grouping populations per locality and suggesting higher divergence between the Corsican pair and the rest of the samples ( Fig. 1e ). Genomic footprints of divergent selection and of association with altitude To identify footprints of divergent selection, including patterns of association with altitude, three distinct BayPass analyses were carried from pool-seq data while taking into account population structure ( Fig. 1f , Fig. 2 ). For each analysis, local scores were computed to identify candidate regions genome-wide by combining SNP-specific statistics at 4,635,467 SNPs with MAF ≥ 10%. Download figure Open in new tab Figure 2. Results of genome-wide scans for signatures of selection at 4,635,467 SNPs (MAF ≥ 10%) genotyped in pool-seq libraries along the 50 chromosomes of the pine processionary genome. For each panel, SNPs located within significant local score regions are depicted in red, with region number and sizes provided in the upper right corner. a. Local scores of the X t X statistic detecting positive selection signatures, irrespective of the elevation levels (ξ = 3, significance level 𝛼 = 0.001). b. Local scores of the C 2 statistic, contrasting low and high altitude groups of populations (ξ = 2, 𝛼 = 0.01). c. Local scores of Bayes Factors ( BF ) measuring the strength of evidence for the linear correlation between population allele frequencies and altitude (ξ = 2, 𝛼 = 0.01). The first analysis focusing on signatures of positive selection ( X t X statistic) identified 721 genomic regions scattered across all 50 chromosomes, with a mean region size of 2.53 kbp (min: 23 bp, max: 134 kbp) and a cumulated size of 1.82 Mbp (ca. 3‰ of the genome, Fig. 2a , Suppl. Table 2). This set of regions overlapped with 227 unique gene models (Suppl. Table 3). The second analysis contrasting high- and low-elevation populations ( C 2 statistic) identified 23 genomic regions scattered across 14 chromosomes, with a mean region size of 3.36 kbp (min: 420 bp, max: 11.4 kbp) and a cumulated size of 77.3 kbp ( Fig. 2b , Suppl. Table 4). These regions overlapped with 11 gene models (Suppl. Table 5). Finally, the third analysis aimed at identifying associations between population allele frequencies and altitude through the computation of Bayes Factors ( BF ), and identified a single 7.43 kbp genomic interval located on chromosome 4 ( Fig. 2c , Suppl. Table 6), which overlapped with two gene models (Suppl. Table 7). No significant gene ontology enrichment was detected in any of the three candidate gene sets. While there was only little overlap between C 2 and X t X regions (13 genomic intervals), all three statistics detected the same major signal between 12 and 13 Mbp on chromosome 4 (Suppl. Fig. 2). A strong putative selective sweep on chromosome 4 To further dissect the signature of selection detected on chromosome 4, site-specific and multilocus windowed pairwise F ST were computed within each locality to identify which population pairs were driving the signal. While little to no divergence was observed in the Spanish and Italian population pairs, divergence was accrued in the Mont Ventoux and Corsican pairs ( Fig. 3 ). Interestingly, for these two pairs, this region was the most differentiated genome-wide between low and high altitude populations, albeit at different levels ( F ST estimates of 0.52 and 0.15 for Corsica and Mont Ventoux, respectively, Suppl. Fig. 3). Parallel divergence signatures between Corsican and Mont Ventoux population pairs could be caused by the introgression of an advantageous haplotype from one locality to the other. To test for this hypothesis, we computed the normalized F 4 statistic (Patterson’s D , Green et al. 2010) along the genome using the poolfstat R package, testing all possible configurations between both Corsican and Mont Ventoux populations, Spanish populations and recently released North African T. pityocampa data as outgroups (Muller et al. 2025). Unfortunately, none of the tested configurations met the treeness criterion genome-wide, suggesting gene flow between populations and limiting our ability to detect putative adaptive introgression with the data at hand (data not shown). We thus focused on the Corsican pair to better characterize patterns of nucleotide diversity and divergence at this locus ( Fig. 4 , see Suppl. Fig. 4, 5 and 6 for results in the other localities). Download figure Open in new tab Figure 3. Pairwise F ST divergence patterns within each locality in the main outlier genomic region on chromosome 4. Plain lines represent 100-SNP sliding window multilocus estimates. Purple blocks indicate candidate X t X regions identified with the local score approach. Dotted lines give the 95 th quantile of the genome-wide F ST distribution for each population pair. a. Spain. b. Italy. c. Mont Ventoux. d. Corsica. Download figure Open in new tab Figure 4. Zoom on the selective sweep signal detected between the two Corsican populations on chromosome 4. a. SNP-specific pairwise F ST values between low and high altitude pools (green: multilocus F ST values computed in sliding 100-SNP windows). Purple blocks indicate candidate X t X regions identified with the local score approach. b. Nucleotide diversity within (π) and absolute divergence between ( d XY ) populations, computed from ind-seq data in 5 kbp non-overlapping windows. c. Tajima’s D computed from pools in 5 kbp non-overlapping windows, with genome-wide means indicated as plain lines. d. Derived allele frequencies in high and e. low altitude pools, with SNPs colored according to their estimated functional impact predicted with SNPeff (yellow: synonymous, orange: non-synonymous, red: putative loss-of-function). f. Standardized sequencing depths computed from ind-seq libraries in 10 kbp windows. Dotted lines at 0.5 and 2 indicate halved and doubled average depths, respectively. g. Results of a screening for structural variation ≥ 500 base pairs using all Corsican ind-seq libraries jointly (DEL: deletion, INV: inversion, DUP: duplication). h. Gene annotation track (cytochrome P450 genes indicated in blue with asterisks). Pairwise F ST estimates between low and high Corsican pools were the highest between 12.2 and 12.8 Mbp, with 69 SNPs displaying F ST values exceeding 0.9, including 14 SNPs alternatively fixed between low and high populations in the pool-seq data ( Fig. 4a ). Close to this region of elevated relative divergence, nucleotide diversity π and absolute divergence d XY computed from ind-seq data were low in both populations between 11.375 and 12.13 Mbp, while flanking regions showed higher absolute divergence levels ( Fig. 4b ). In this same low nucleotide diversity region, Tajima’s D values were negative for both populations ( Fig. 4c ), which both showed deficits of SNPs at intermediate frequencies ( Fig. 4d , Fig. 4e ). Altogether, these observations (low nucleotide diversity and negative Tajima’s D ) are in agreement with a scenario of selective sweep, and the functions of genes in the region may provide insights into the factors that have driven this putative selective event. Functional insights from the sweep region highlight detoxification mechanisms Narrowing the focus to the 12-13 Mbp interval on chromosome 4, which encompasses the candidate X t X regions, we recovered 20 gene models, including a cluster of 5 cytochrome P450 monooxygenase genes from the CYP6B gene family (between 12.08 and 12.32 Mbp, with a sixth CYP6B gene outside the interval around 11.33 Mbp, Fig. 4 and Suppl. Fig. 7). Importantly, the expert annotation of detoxification genes in the pine processionary genome confirmed all P450 gene models predicted within this region were complete and showed alternatively spliced forms ( Gautier et al. 2025 ), with expression support from RNA-seq data across multiple populations and developmental stages (Suppl. Fig. 8). Both individual sequencing depth profiles and short-read-based structural variation screening suggest a complex structural landscape in the vicinity of the cytochrome P450 cluster ( Fig. 4f and 4g). Despite our stringent SNP filtering scheme, such indirect evidence for structural variation in the genomic region makes the fine mapping of putative causal substitutions within P450 genes hazardous. The 12-13 Mbp region also contains the voltage-gated sodium channel gene para , which is located 23 kbp away from the P450 cluster, in-between two candidate X t X regions (ca. 12.34 Mbp, Fig. 4 and Suppl. Fig. 7). Of note, no evidence of adaptive differentiation or altitude-associated variation was detected for SNPs within para (Suppl. Fig. 7, synonymous variants in yellow). Among other genes in this 1 Mbp interval, an endonuclease-reverse transcriptase located only 3 kbp upstream of the P450 cluster showed strong signals of adaptive divergence (including one highly differentiated non-synonymous substitution, log10( PXtX ) = 12.0, F Corsica = 0.92). However, a high number of non-synonymous substitutions (Suppl. Fig. 7) together with very low expression levels across multiple populations and developmental stages (Suppl. Fig. 8) suggest that this gene is undergoing pseudogenization. Finally, a group of four genes located between 12.6 and 12.8 Mbp also displayed moderate X t X signatures and strong signals of differentiation between the Corsican population pairs. These four genes included two oxidoreductase enzymes (LOC113503460), along with homologs of Rabconnectin-3A and Tempura , both of which are involved in the Notch signaling pathway in Drosophila ( Tuttle et al. 2014 ; Charng et al. 2014 ). Discussion A strong selective sweep signature suggestive of insecticide resistance Our replicated genome-wide scan across four altitudinal gradients in the pine processionary moth uncovered a strong and repeated footprint of selection across a 1.5 Mb region located on chromosome 4. With a marked reduction in nucleotide diversity, negative Tajima’s D , and near-fixation of alternative alleles between Corsican populations, this region displays the classical genomic hallmarks of a recent selective sweep. Remarkably, the putative sweep encompasses a cluster of CYP6B cytochrome P450 monooxygenase genes and the voltage-gated sodium channel gene para , both widely implicated in insecticide detoxification and target-site resistance in insects ( Puinean et al. 2010 ; Liu 2015 ; Kaduskar et al. 2022 ; Ngegba et al. 2025 ; Nelson et al. 1996 ; Pu & Chung 2024 ). The recurrence of this selective footprint in two geographically distant population pairs (Corsica and Mont Ventoux, albeit less clearly, Suppl. Fig. 5), together with the lack of evidence for recent introgression between these localities and the strong geographic structure detected among localities, suggest independent adaptive events rather than the spread of a single advantageous haplotype. Recurrent, independent sweeps at the same locus have been documented in several insect pests evolving resistance through gene amplification, regulatory changes, or target-site substitutions ( Bass & Field 2011 ; Pélissié et al. 2022 ; Chen et al. 2025 ). While there is little support in our dataset for point mutations driving reduced sensitivity to insecticides in the para gene, coverage profiles and short-read structural-variant scans reveal a complex structural landscape surrounding the CYP6B cluster, potentially including duplications and/or rearrangements. Copy number expansion of P450 genes is a well-characterized mechanism leading to xenobiotic resistance, where additional copies of detoxification genes increase transcriptional output or provide substrates for neofunctionalization ( Scott 1999 ; Bass & Field 2011 ). Moreover, recent mechanistic studies in Spodoptera exigua showed that both cis- and trans-regulatory changes synergistically increase P450 expression to confer resistance to both organophosphate and pyrethroid insecticides ( Hu et al. 2025 ), suggesting that both structural and regulatory evolution could drive the putative selective sweep detected in the vicinity of the CYP6B cluster. The role of CYP6 genes in pesticide detoxification and resistance is well established (Nauen et al. 2022). In Lepidoptera, their overexpression is classically associated with resistance to synthetic pyrethroids (e.g., Bautista et al. 2009; do Nascimento et al. 2023; Xu et al. 2020; Perini et al. 2021). Pyrethroids act as neurotoxic agents that disrupt the insect’s nervous system, mostly targeting voltage-gated sodium channels in insect nerve cells, causing paralysis and death. In France, the use of pyrethroids compounds (e.g., deltamethrin, bifenthrin, cypermethrin, permethrin, fenvalerate/esfenvalerate, lambda-cyhalothrin) was authorized during several decades in forests notably for controlling the pine processionary moth, but their use has progressively been banned during the 2000s and 2010s in France and in Europe. Pine processionary moth management now relies notably on Bacillus thuringiensis extracts ( Bt toxin) at large scales or on caterpillar traps at local scales. As CYP6 genes have so far not been shown to be implicated in resistance to Bt toxin (Heckel et al. 2007; Deqin Hu et al. 2025 ), the putative selective sweep identified likely predates the adoption of Bt -based control and could have been caused by the previous use of pyrethroids. The maintenance of resistance alleles in natural populations following the ban on pyrethroids may be explained primarily by the relatively recent effective ban of these compounds and by a potentially low cost of the resistance mechanism. Interestingly, Calla et al. (2021) investigated insecticide resistance in another moth, the navel orangeworm Amyelois transitella , which is a major almond pest in US orchards. Comparing populations with different levels of resistance to bifenthrin, a pyrethroid, they also detected a putative hard selective sweep covering 1.3 Mbp and encompassing the same candidate P450 cluster and para gene, suggesting gene order conservation between A. transitella and T. pityocampa despite >100 MY of divergence (the syntenic block extends beyond kruppel and ariadne , shown on Fig. 4h ). In the future, combining long read sequencing and gene expression profiling across populations with differing exposure histories and/or controlled experiments under different insecticide exposures could help confirm the presence of structural variants and study their phenotypic effects among pine processionary populations. Tree secondary metabolites as potential drivers of detoxification genes Another potential selective pressure driving P450 expansion is resistance to natural tree defense compounds such as α-pinene, a monoterpene present in pine tissues, including the needles that caterpillars feed on (Kubeczka & Schultze 1987). P450s (notably CYP4 and CYP6 gene families) are central to pinene detoxification in Coleoptera such as bark beetles ( Dai et al. 2014 , 2015 ; Chiu et al. 2019 ) and longhorn beetles ( He et al. 2022 ). In these species, the expression of detoxification genes in response to α-pinene and other monoterpenes has been shown to vary among populations exposed to different host trees or terpene concentrations ( He et al. 2022 ; Torres-Banda et al. 2022 ). While evidence is less abundant in Lepidoptera, α-pinene has been shown to induce overexpression of CYP6B2 in the cotton bollworm Helicoverpa armigera ( Ranasinghe et al. 1997 ), and exposure to α-pinene led to both oxidative stress and nervous system disruption in the flour moth Ephestia kuehniella ( Shahriari et al. 2018 ), increasing the enzymatic activity of superoxide dismutases and glutathione S-transferases (P450s were not assayed in the study). In light of these observations, the selection footprint we detected could reflect response to natural α-pinene exposure. Under this hypothesis, selection may be acting via differences in α-pinene concentrations across pine trees along the altitudinal gradient, maybe related to variations in pine species assemblages (the α-pinene to β-pinene ratio is variable among pine species, Malhocká et al. 2024 ). Both insecticide resistance and terpene detoxification could act concurrently, exerting selection in the same direction, on the same detoxification genes and/or within the same P450 cluster. Analyses of altitudinal variation in α-pinene concentrations across study sites, coupled with transcriptional profiling of populations exposed to different tree defense compounds, would further test this hypothesis. Weak signatures of altitudinal adaptation in the pine processionary moth The genomics of altitudinal adaptation are increasingly studied in Lepidoptera, revealing both parallel and polygenic genomic architectures. At a larger evolutionary timescale compared to our study, genome-wide scans in Heliconius erato and H. melpomene along replicated altitudinal transects identified multiple divergent regions associated with altitude, showing that adaptation often draws on standing variation or adaptive introgression from related taxa and involves moderate-effect loci distributed across the genome ( Montejo-Kovacevich et al. 2022 ). In the same system, a related genome-wide association study of wing aspect ratio showed that altitude-related morphological shifts are heritable and highly polygenic ( Montejo-Kovacevich et al. 2019 ). In our study, outside the major signal on chromosome 4, we found only a few genomic regions displaying elevated levels of divergence between low- and high-altitude groups of populations and/or association with altitude ( C 2 and BF statistics, respectively), for which gene content analysis could not provide any clear link between gene function and selective pressures possibly experienced in the wild (see Suppl. Tables 3, 5 and 7). Of note, while the X t X analysis aimed at detecting local adaptation signatures uncovered 721 genomic regions totaling 1.82 Mbp, the C 2 analysis only identified 23 genomic regions totaling 77.3 kbp, thus suggesting that altitude may not exert a major selection pressure, at least within our study area. Alternatively, assuming altitude is a major stressor for the pine processionary moth, such sparse signals can emerge under several non-mutually exclusive scenarios, including polygenic adaptation governed by many small-effect loci ( Pritchard et al. 2010 ), substantial phenotypic plasticity buffering selection ( Forsman 2015 ), heterogeneous responses to selection across localities ( Barghi et al. 2020 ), or interference from stronger, unrelated selective forces ( McCulloch & Waters 2023 ). Adaptation to altitude is also influenced by a combination of factors, including pine processionary phenology ( Martin et al. 2022 ), temperature, air pressure differences, and changes in habitat structure and composition ( Sømme 1989 ). Adding to this complexity, organismal responses to these factors can be correlated, and disentangling their individual effects from our observational data alone is challenging. While the methods we used to recover divergence and association signals from genome-wide polymorphism data show modest to good performance to detect polygenic adaptation in simulations ( Gautier 2015 ), pinpointing small-effect loci remains a difficult task ( Barghi et al. 2020 ). Additionally, altitudes were not fully balanced between localities in our sampling scheme ( Table 1 ), owing to difficulties in localizing and/or sampling populations and/or nests, and such imbalance could hamper the power of our association analyses (e.g., the elevation of the low-elevation Spanish population is higher than that of the high-elevation Corsican population). Tangled selective pressures in natural populations Our study exemplifies the difficulty of teasing apart ecological and anthropogenic selection in natural populations. In many arthropods, insecticide resistance evolves orders of magnitude faster than climatic adaptation because pesticide exposure imposes extreme and episodic selective pressures ( McCulloch & Waters 2023 ). Strong selective sweeps at detoxification or target-site loci can therefore interfere with and mask subtler polygenic responses to temperature or altitude. In addition, linkage and pleiotropy may blur these processes: selection on detoxification genes could indirectly alter physiological traits related to thermal or oxidative stress, while climate-driven selection could modify resistance phenotypes (e.g., Rostant et al. 2017 ). Spatial heterogeneity exposure further complicates inference. Localized pesticide use can generate parallel sweeps in distant populations facing similar treatment regimes but leave other populations unaffected, creating an apparent mismatch between ecological and genomic signals. Disentangling these effects requires a wealth of environmental data, along with controlled assays, which are not realistic in most settings. However, this work also shows that careful interpretation of genomic data can point towards unexpected selection pressures. From a management perspective, the potential for resistance evolution in T. pityocampa , which is in theory not subject to routine treatment, highlights the need for caution and control in chemical interventions carried in natural environments. Methods Study populations and samples Sampling was performed between October 2011 and March 2014 for four pairs of populations along altitudinal gradients located in the Spanish Sierra Nevada, the Italian Alps, Corsica and Mont Ventoux (France; Table 1 and Figure 1 ). For each of the four pairs, individuals were sampled for both an upper slope population and a lower slope population (altitudes are given in Table 1 ). For Corsican and Italian populations, larvae were directly collected from the nests during winters between 2011 and 2014 ( Table 1 ). For Spain and Mont Ventoux, male adults were sampled using pheromone traps during summers 2012 and 2013 ( Table 1 ). All individuals were preserved in 95% ethanol after sampling. To decrease relatedness among individuals both within populations and within pools (see below), a single caterpillar per colony and per tree was collected for both Italian and Corsican populations, while for pheromone traps, individuals were collected using several traps deployed over most of the emergence period (six to eight weeks). Lab work and whole-genome sequencing DNA from each individual was extracted using Qiagen’s commercial DNeasy Blood & Tissue Kit. Population samples ranged from 26 to 40 individuals per population ( Table 1 ). For pooled sequencing, DNA extractions were mixed in equimolar proportions for each population. For individual sequencing, the DNA of 8 to 16 individuals per population were selected, for a total of 89 individuals. A total of 97 libraries (89 ind-seq + eight pool-seq) were constructed using Illumina’s TruSeq nano DNA kit following manufacturer instructions. Paired-end sequencing was carried out at the Montpellier GenomiX (MGX) platform on Illumina HiSeq 2500 (pool-seq samples, read length: 125 nucleotides) and Illumina NovaSeq 6000 (ind-seq samples, read length: 150 nucleotides), targeting mean read depths of 60× and 20× for pool-seq and ind-seq libraries, respectively. Bioinformatics Unless stated otherwise, all analyses and data visualization were carried in R v4.5.1 ( R Core Team 2019 ). Data quality was checked with FastQC v0.12.0 ( Andrews 2010 ) before and after trimming of raw reads using fastp v0.23.2 ( Chen et al. 2018 ). Trimmed reads were aligned against the chromosome-level pine processionary genome assembly Tpit_2.1 ( Gautier et al. 2025 ) using bwa mem v0.7.18 ( Li 2013 ) and duplicates were removed with sambamba v1.0.1 ( Tarasov et al. 2015 ). Insert sizes were collected from BAM files using Picard v3.0.0 ( http://broadinstitute.github.io/picard ), and sequencing depths were tallied for each library in 10 kbp windows along the genome with mosdepth v0.3.6 ( Pedersen & Quinlan 2018 ). SNP calling was performed jointly on all libraries using the haplotype-aware caller FreeBayes v1.3.6 ( Garrison & Marth 2012 ) with the following parameters: --pooled-continuous -m 20 -q 20 -p 2 -n 3 --min-coverage 10 --limit-coverage 300 -E - 1 --max-complex-gap -1 --haplotype-length -1 . Multi-nucleotide variants were decomposed with vt v0.57721 ( Tan et al. 2015 ), and only biallelic sites displaying qualities ≥ 30 and supported by at least one forward and one reverse reads across all ind-seq and pool-seq libraries were retained using BCFtools v1.17 ( Danecek et al. 2021 ). The resulting VCF file was split into pool-seq ( n = 8 samples) and ind-seq VCF files ( n = 89 samples) to apply library-specific filtering schemes. For the pool-seq VCF file, sites were filtered using the poolfstat v3.0.0 R package ( Gautier et al. 2022 ) to retain SNPs displaying overall minor allele frequencies (MAF) computed on read counts ≥ 0.05, at least four observations of the minor allele across all pools, sequencing depths comprised between 10× and 200× in each pool, and discarding sites located at less than five base pairs from an indel. For the ind-seq VCF file, since samples contained larvae for which sex cannot be determined based on morphology, Z-chromosome-to-autosomal sequencing depth ratios were computed for each ind-seq sample. Females are heterogametic in Lepidoptera and should display sequencing depth ratios around 0.5 (one Z copy for two autosomal copies), compared to unity for males. The depth ratio analysis identified 66 males and 23 females from the ind-seq data, with all individuals sampled using pheromone traps assigned as males, as expected ( Table 1 , Suppl. Table 1, Suppl. Fig. 1). After this bioinformatic sex determination step, biallelic sites displaying depths below 8× and above twice the mean library depth for each individual were coded as missing using BCFtools. The lower filter was adjusted for Z-linked sites in females (allowing minimum 4×). Sites displaying more than 25% missing data across all individuals were discarded, along with sites which did not fulfil the global MAF ≥ 0.05 criterion (computed across all ind-seq libraries). Additionally, in order to compute unbiased nucleotide diversity and divergence metrics, a VCF file containing both variant and invariant sites (hereafter, all sites) was generated solely from ind-seq data using BCFtools, following recommendations from Korunes & Samuk (2021) . Population structure, divergence and nucleotide diversity metrics Population structure was first documented from ind-seq data using a principal component analysis (PCA) run on 84,443 linkage-pruned autosomal genotypes with plink v2.0 (pruning carried with parameters --indep-pairwise 250kb 1 0.2, Chang et al. 2015 ). Additionally, unsupervised clustering was carried on the same dataset with Admixture v1.3.0 ( Alexander et al. 2009 ) using default parameters and ranging ancestry components K from 2 to 12, with 10 replicates per K value. Results were aggregated with StructureSelector ( Li & Liu 2018 ) and the most likely number of ancestry components was estimated following Puechmaille (2016) , owing to unbalanced sample sizes across populations. Nucleotide diversity (π) and pairwise divergence metrics ( F ST and d XY) were computed from the all sites ind-seq VCF file in 5 kbp windows along the genome using pixy v1.2.11 ( Korunes & Samuk 2021 ). Identification of candidate genomic regions associated with altitude Genomic regions displaying altitudinal variations in allele frequencies were identified using the hierarchical Bayesian framework implemented in BayPass v3.0.0 ( Gautier 2015 ), which explicitly accounts for the covariance structure among population allele frequencies resulting from the shared history of samples. Three distinct analyses were carried out for each SNP of the pool-seq dataset under BayPass’ standard covariate model, specifying pool haploid sample sizes (assuming all individuals within a pool are male and carry two copies of the Z chromosome, which is conservative) and using default options for the importance sampling approximation. First, the X t X statistic ( Günther & Coop 2013 ) was computed to detect signals of positive selection in any population, irrespective of elevation levels. Second, the C 2 statistic ( Olazcuaga et al. 2020 ) was computed to contrast allele frequencies between the low- and high-altitude groups of populations (i.e., altitude is treated as a binary trait, with two groups of four populations each). Third, Bayes factors ( BF ) were computed to measure the strength of evidence for the linear association between population allele frequencies and altitude in meters (scaled by default in BayPass, so that 𝜇 = 0 and σ 2 = 1). For all three analyses, autosomal and Z-linked SNPs were treated separately, and to alleviate the computational burden, the full dataset was thinned into subsets containing ca. 75,000 SNPs each (leading to 185 and 6 subsets for autosomal and Z compartments, respectively). The local score approach developed by Fariello et al. (2017) was then applied to the statistics computed with BayPass to account for linkage. This approach implemented in BayPass’ companion R script combines statistics at neighboring SNPs with global MAF ≥ 10% along each chromosome to detect significant regions based on the user-defined significance threshold ξ. To better prioritize candidate genomic regions, distinct thresholds were applied to the different statistics: ξ = 2 for both C 2 and BF (significance level 𝛼 = 0.01), versus ξ = 3 for X t X (𝛼 = 0.001). For each analysis, gene models located within candidate genomic regions were extracted from the genome annotation file and gene ontology enrichments were tested using the topGO R package v2.60.1 ( Alexa & Rahnenfuhrer 2025 ). Fine-scale characterization of the main candidate selective sweep region Genome-wide analyses using BayPass recovered a strong signal on chromosome 4, which was further characterized at a finer scale. To identify which population pair(s) drove the signal in the genomic region, pairwise F ST metrics were computed between low and high populations within each locality from pool-seq data. Then, focusing on the Corsican population pair, which showed the highest relative divergence in the region (see Results), and in addition to π and d XY metrics previously computed from ind-seq data, Tajima’s D ( D T) was computed from pool-seq data in 5 kbp windows using npstat v1 with default options ( Ferretti et al. 2013 ). All SNPs passing previous filters in the pool-seq VCF file were polarized by setting as ancestral the major allele from the two Spanish populations, as Spain is assumed to be a major glacial refugium of Western European pine processionary populations ( Kerdelhué et al. 2009 ). These polymorphisms were further annotated with SNPeff v4.0 ( Cingolani 2022 ) using the Tpit_2.1 gff annotation file ( Gautier et al. 2025 ). Finally, to better characterize the structural landscape in the candidate genomic region, sequencing depths were tallied for each Corsican ind-seq library in 10 kbp windows using mosdepth v0.3.6 ( Pedersen & Quinlan 2018 ). Additionally, evidence of structural variation (SV) was gathered from the same dataset in the region using smoove v0.2.8 ( https://github.com/brentp/smoove ), a wrapper which implements internal read filtering steps and combines the Lumpy SV caller ( Layer et al. 2014 ) to Duphold for depth-based SV annotation ( Pedersen & Quinlan 2019 ). Detecting structural variation from short-read data is a complex endeavor (see, e.g., Olson et al. 2023 ), and results from this last analysis should thus be interpreted cautiously. Data accessibility The data underlying this article are deposited in the European Nucleotide Archive and will be made available pending acceptance of the manuscript. Sampling locations are provided in Table 1 and sample metadata are provided in Suppl. Table 1. Acknowledgments This work was funded through both the ANR ANR-22-CE02-0009 LOADEXP and INRAE exploratory projects to C.P. and ANR-10-JCJC-1705-01 GENOPHENO to C.K. We are grateful to Aurélie Ducasse, Renaud Vitalis, Vincent Danneville, Jean-Michel Gigleux, Pikria Doinjashvili and Cécile Robin for their assistance with project management. We acknowledge our colleagues Andrea Battisti (Univ. Padova, Italy), José A. Hódar and Rut Aspizua (Univ. Granada, Spain), Jérôme Rousselet, Alain Roques and Francis Goussard (URZF, INRAE Orléans Val de Loire, France) and Jean-Claude Martin, Anne-Sophie Brinquin, Marianne Correard, René Mazet and Estelle Morel (UEFM, INRAE PACA, France) for providing samples from Italy, Spain, Corsica and Mont Ventoux, respectively. Sequencing services were provided by the MGX facility ( https://mgx.cnrs.fr ). We acknowledge the GenSeq platform (Université Montpellier) for granting us access to their facilities for library preparation. We are grateful to the genotoul bioinformatics platform Toulouse Occitanie (Bioinfo Genotoul, https://doi.org/10.15454/1.5572369328961167E12 ) and to the CBGP cluster (Nathalie Vieira and Sylvain Piry) for providing computing and storage resources. Funder Information Declared Agence Nationale de la Recherche , ANR-22-CE02-0009 , ANR-10-JCJC-1705-01 INRAE , Recurrent funding Bibliography ↵ Alexa A , Rahnenfuhrer J. 2025 . topGO: Enrichment Analysis for Gene Ontology . https://bioconductor.org/packages/topGO (Accessed August 19, 2025). ↵ Alexander DH , Novembre J , Lange K . 2009 . Fast model-based estimation of ancestry in unrelated individuals . Genome Res . 19 : 1655 – 1664 . OpenUrl Abstract / FREE Full Text ↵ Andrews S. 2010 . FastQC: a quality control tool for high throughput sequence data . https://www.bioinformatics.babraham.ac.uk/projects/fastqc . ↵ Barghi N , Hermisson J , Schlötterer C . 2020 . Polygenic adaptation: a unifying framework to understand positive selection . Nat. Rev. Genet . 21 : 769 – 781 . OpenUrl CrossRef PubMed ↵ Bass C , Field LM . 2011 . Gene amplification and insecticide resistance . Pest Manag. Sci . 67 : 886 – 890 . OpenUrl CrossRef PubMed Web of Science ↵ Battisti A et al. 2005 . Expansion of geographic range in the pine processionary moth caused by increased winter temperatures . Ecol. Appl . 15 : 2084 – 2096 . OpenUrl ↵ Bellard C , Bertelsmeier C , Leadley P , Thuiller W , Courchamp F . 2012 . Impacts of climate change on the future of biodiversity: Biodiversity and climate change . Ecol. Lett . 15 : 365 – 377 . OpenUrl CrossRef PubMed ↵ Bernatchez L , Ferchaud A-L , Berger CS , Venney CJ , Xuereb A . 2024 . Genomics for monitoring and understanding species responses to global climate change . Nat. Rev. Genet . 25 : 165 – 183 . OpenUrl CrossRef PubMed ↵ Bonhomme M et al. 2010 . Detecting selection in population trees: The Lewontin and Krakauer test extended . Genetics . 186 : 241 – 262 . OpenUrl Abstract / FREE Full Text ↵ Bourgeois YXC , Warren BH . 2021 . An overview of current population genomics methods for the analysis of whole-genome resequencing data in eukaryotes . Mol. Ecol . 30 : 6036 – 6071 . OpenUrl CrossRef ↵ Bourougaaoui A , Ben Jamâa ML , Robinet C . 2021 . Has North Africa turned too warm for a Mediterranean forest pest because of climate change? Clim. Change . 165 : 1 – 20 . OpenUrl ↵ Calla B et al. 2021 . Selective sweeps in a nutshell: The genomic footprint of rapid insecticide resistance evolution in the almond agroecosystem . Genome Biol. Evol . 13 : evaa234 . OpenUrl PubMed ↵ Chang CC et al. 2015 . Second-generation PLINK: rising to the challenge of larger and richer datasets . Gigascience . 4 : 7 . OpenUrl CrossRef PubMed ↵ Charng W-L et al. 2014 . Drosophila Tempura, a novel protein prenyltransferase α subunit, regulates notch signaling via Rab1 and Rab11 . PLoS Biol . 12 : e1001777 . OpenUrl CrossRef PubMed ↵ Chen I-C , Hill JK , Ohlemüller R , Roy DB , Thomas CD . 2011 . Rapid range shifts of species associated with high levels of climate warming . Science . 333 : 1024 – 1026 . OpenUrl Abstract / FREE Full Text ↵ Chen L et al. 2025 . Population genomics reveal multiple independent origins of pesticide resistance in the polyphagous pest , Tetranychus urticae. Commun. Biol . 8 : 1197 . OpenUrl PubMed ↵ Chen S , Zhou Y , Chen Y , Gu J . 2018 . fastp: an ultra-fast all-in-one FASTQ preprocessor . Bioinformatics . 34 : i884 – i890 . OpenUrl CrossRef PubMed ↵ Chiu CC , Keeling CI , Bohlmann J . 2019 . The cytochrome P450 CYP6DE1 catalyzes the conversion of α-pinene into the mountain pine beetle aggregation pheromone trans-verbenol . Sci. Rep . 9 : 1477 . OpenUrl CrossRef PubMed ↵ Cingolani P . 2022 . Variant annotation and functional prediction: SnpEff . Methods Mol. Biol . 2493 : 289 – 314 . OpenUrl CrossRef PubMed ↵ Dai L et al. 2015 . Cytochrome P450s from the Chinese white pine beetle, Dendroctonus armandi (Curculionidae: Scolytinae): Expression profiles of different stages and responses to host allelochemicals . Insect Biochem. Mol. Biol . 65 : 35 – 46 . OpenUrl PubMed ↵ Dai L et al. 2014 . Two CYP4 genes of the Chinese white pine beetle, Dendroctonus armandi (Curculionidae: Scolytinae), and their transcript levels under different development stages and treatments: Dendroctonus armandi . Insect Mol. Biol . 23 : 598 – 610 . OpenUrl PubMed ↵ Danecek P et al. 2021 . Twelve years of SAMtools and BCFtools . Gigascience . 10 . doi: 10.1093/gigascience/giab008 . OpenUrl CrossRef PubMed ↵ Dudaniec RY , Yong CJ , Lancaster LT , Svensson EI , Hansson B . 2018 . Signatures of local adaptation along environmental gradients in a range-expanding damselfly (Ischnura elegans) . Mol. Ecol . 27 : 2576 – 2593 . OpenUrl CrossRef ↵ Fariello MI et al. 2017 . Accounting for linkage disequilibrium in genome scans for selection without individual genotypes: The local score approach . Mol. Ecol . 26 : 3700 – 3714 . OpenUrl CrossRef ↵ Ferretti L , Ramos-Onsins SE , Pérez-Enciso M . 2013 . Population genomics from pool sequencing . Mol. Ecol . 22 : 5561 – 5576 . OpenUrl CrossRef PubMed Web of Science ↵ Forester BR , Lasky JR , Wagner HH , Urban DL . 2018 . Comparing methods for detecting multilocus adaptation with multivariate genotype-environment associations . Mol. Ecol . 27 : 2215 – 2233 . OpenUrl CrossRef ↵ Forsman A . 2015 . Rethinking phenotypic plasticity and its consequences for individuals, populations and species . Heredity . 115 : 276 – 284 . OpenUrl CrossRef PubMed ↵ Futschik A , Schlötterer C . 2010 . The next generation of molecular markers from massively parallel sequencing of pooled DNA samples . Genetics . 186 : 207 – 218 . OpenUrl Abstract / FREE Full Text ↵ Garrison E , Marth G. 2012 . Haplotype-based variant detection from short-read sequencing . arXiv [q-bio.GN]. http://arxiv.org/abs/1207.3907 . ↵ Gautier M et al. 2025 . A chromosome-level assembly of the pine processionary moth ( Thaumetopoea pityocampa ) genome . bioRxiv . 2025.06.06.658235. doi: 10.1101/2025.06.06.658235 . OpenUrl Abstract / FREE Full Text ↵ Gautier M . 2015 . Genome-Wide Scan for Adaptive Divergence and Association with Population-Specific Covariates . Genetics . 201 : 1555 – 1579 . OpenUrl Abstract / FREE Full Text ↵ Gautier M , Vitalis R , Flori L , Estoup A . 2022 . ƒ-statistics estimation and admixture graph construction with Pool-Seq or allele count data using the R package poolfstat . Mol. Ecol. Resour . 22 : 1394 – 1416 . OpenUrl CrossRef PubMed ↵ Günther T , Coop G . 2013 . Robust identification of local adaptation from allele frequencies . Genetics . 195 : 205 – 220 . OpenUrl Abstract / FREE Full Text ↵ He X et al. 2022 . Transcriptome profiling and RNA interference reveals relevant detoxification genes in Monochamus alternatus response to (+)-α-pinene . J. Appl. Entomol . 146 : 823 – 837 . OpenUrl ↵ Hill JK , Griffiths HM , Thomas CD . 2011 . Climate change and evolutionary adaptations at species’ range margins . Annu. Rev. Entomol . 56 : 143 – 159 . OpenUrl CrossRef PubMed Web of Science ↵ Hu B et al. 2025 . Two independent regulatory mechanisms synergistically contribute to P450-mediated insecticide resistance in a lepidopteran pest, Spodoptera exigua . BMC Biol . 23 : 122 . OpenUrl PubMed ↵ Inouye DW . 2022 . Climate change and phenology . Wiley Interdiscip. Rev. Clim. Change . 13 . doi: 10.1002/wcc.764 . OpenUrl CrossRef ↵ Kaduskar B et al. 2022 . Reversing insecticide resistance with allelic-drive in Drosophila melanogaster . Nat. Commun . 13 : 291 . OpenUrl CrossRef PubMed ↵ Kerdelhué C et al. 2009 . Quaternary history and contemporary patterns in a currently expanding species . BMC Evol. Biol . 9 : 220 . OpenUrl CrossRef PubMed Kerdelhué C , Rossi J-P , Rousselet J . 2023 . Expansion fronts of the pine processionary moth Thaumetopoea pityocampa from 2005 to 2021 . doi: 10.57745/QG9P9R . OpenUrl CrossRef ↵ Korunes KL , Samuk K . 2021 . pixy: Unbiased estimation of nucleotide diversity and divergence in the presence of missing data . Mol. Ecol. Resour . 21 : 1359 – 1368 . OpenUrl CrossRef PubMed ↵ Layer RM , Chiang C , Quinlan AR , Hall IM . 2014 . LUMPY: a probabilistic framework for structural variant discovery . Genome Biol . 15 : R84 . OpenUrl CrossRef PubMed ↵ Li H. 2013 . Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM . arXiv [q-bio.GN]. http://arxiv.org/abs/1303.3997 (Accessed October 10, 2024). ↵ Liu N . 2015 . Insecticide resistance in mosquitoes: impact, mechanisms, and research directions . Annu. Rev. Entomol . 60 : 537 – 559 . OpenUrl CrossRef PubMed ↵ Li Y-L , Liu J-X . 2018 . StructureSelector: A web-based software to select and visualize the optimal number of clusters using multiple methods . Mol. Ecol. Resour . 18 : 176 – 177 . OpenUrl CrossRef PubMed ↵ Lotterhos KE , Whitlock MC . 2015 . The relative power of genome scans to detect local adaptation depends on sampling design and statistical method . Mol. Ecol . 24 : 1031 – 1046 . OpenUrl CrossRef PubMed ↵ Malhocká A , Švábová M , Havelcová M . 2024 . The monoterpenes as a characteristic marker for pine species distinction: a chemotaxonomic study from the Czech Republic . J. Essent. Oil-Bear. Plants . 27 : 870 – 886 . OpenUrl ↵ Martin J-C , Mesmin X , Buradino M , Rossi J-P , Kerdelhué C . 2022 . Complex drivers of phenology in the pine processionary moth: Lessons from the past . Agric. For. Entomol . 24 : 247 – 259 . OpenUrl CrossRef ↵ McCulloch GA , Waters JM . 2023 . Rapid adaptation in a fast-changing world: Emerging insights from insect genomics . Glob. Chang. Biol . 29 : 943 – 954 . OpenUrl PubMed ↵ Montejo-Kovacevich G et al. 2019 . Altitude and life-history shape the evolution of Heliconius wings . Evolution . 73 : 2436 – 2450 . OpenUrl CrossRef PubMed ↵ Montejo-Kovacevich G et al. 2022 . Repeated genetic adaptation to altitude in two tropical butterflies . Nat. Commun . 13 : 4676 . OpenUrl CrossRef PubMed ↵ Nelson DR et al. 1996 . P450 superfamily: update on new sequences, gene mapping, accession numbers and nomenclature . 6 : 1 – 42 . OpenUrl ↵ Ngegba PM , Khalid MZ , Jiang W , Zhong G . 2025 . An overview of insecticide resistance mechanisms, challenges, and management strategies in Spodoptera frugiperda . Crop Prot . 197 : 107322 . OpenUrl ↵ Olazcuaga L et al. 2020 . A whole-genome scan for association with invasion success in the fruit fly Drosophila suzukii using contrasts of allele frequencies corrected for population structure . Mol. Biol. Evol . 37 : 2369 – 2385 . OpenUrl CrossRef PubMed ↵ Olson ND et al. 2023 . Variant calling and benchmarking in an era of complete human genome sequences . 24 : 464 – 483 . OpenUrl ↵ Parmesan C , Yohe G . 2003 . A globally coherent fingerprint of climate change impacts across natural systems . Nature . 421 : 37 – 42 . OpenUrl CrossRef PubMed Web of Science ↵ Pedersen BS , Quinlan AR . 2019 . Duphold: scalable, depth-based annotation and curation of high-confidence structural variant calls . Gigascience . 8 . doi: 10.1093/gigascience/giz040 . OpenUrl CrossRef PubMed ↵ Pedersen BS , Quinlan AR . 2018 . Mosdepth: quick coverage calculation for genomes and exomes . Bioinformatics . 34 : 867 – 868 . OpenUrl CrossRef PubMed ↵ Pélissié B et al. 2022 . Genome resequencing reveals rapid, repeated evolution in the Colorado potato beetle . Mol. Biol. Evol . 39 : msac016 . OpenUrl CrossRef PubMed ↵ Pritchard JK , Pickrell JK , Coop G . 2010 . The genetics of human adaptation: hard sweeps, soft sweeps, and polygenic adaptation . Curr. Biol . 20 : R208 – 15 . OpenUrl CrossRef PubMed Web of Science ↵ Puechmaille SJ . 2016 . The program structure does not reliably recover the correct population structure when sampling is uneven: subsampling and new estimators alleviate the problem . Mol. Ecol. Resour . 16 : 608 – 627 . OpenUrl CrossRef PubMed ↵ Puinean AM et al. 2010 . Amplification of a cytochrome P450 gene is associated with resistance to neonicotinoid insecticides in the aphid Myzus persicae . PLoS Genet . 6 : e1000999 . OpenUrl CrossRef PubMed ↵ Pu J , Chung H . 2024 . New and emerging mechanisms of insecticide resistance . Curr. Opin. Insect Sci . 63 : 101184 . OpenUrl PubMed ↵ Ranasinghe C , Headlam M , Hobbs AA . 1997 . Induction of the mRNA for CYP6B2, a pyrethroid inducible cytochrome p450, inHelicoverpa armigera (Hubner) by dietary monoterpenes . Arch. Insect Biochem. Physiol . 34 : 99 – 109 . OpenUrl CrossRef ↵ R Core Team . 2019 . R: A Language and Environment for Statistical Computing . Vienna, Austria https://www.R-project.org/ . ↵ Rivière J. 2011 . Les chenilles processionnaires du pin : évaluation des enjeux de santé animale . École vétérinaire d’Alfort https://theses.vet-alfort.fr/telecharger.php?id=1356 (Accessed August 20, 2025). ↵ Robinet C et al. 2025 . Déplacement de la processionnaire du pin, Thaumetopoea pityocampa . In: Invasion et expansion d’insectes bioagresseurs forestiers: Quels risques pour la forêt française dans le contexte des changements globauxRobinet, C, Saintonge, F-X, Tassus, X, & Brault, S, editors. Éditions Quae pp. 62 – 73 . ↵ Robinet C , Laparie M , Rousselet J . 2015 . Looking beyond the large scale effects of global change: Local phenologies can result in critical heterogeneity in the pine processionary moth . Front. Physiol . 6 : 334 . OpenUrl PubMed ↵ Rossi J-P et al. 2025 . Warmer and brighter winters than before: Ecological and public health challenges from the expansion of the pine processionary moth (Thaumetopoea pityocampa) . Sci. Total Environ . 978 : 179470 . OpenUrl CrossRef PubMed ↵ Rostant WG et al. 2017 . Pleiotropic effects of DDT resistance on male size and behaviour . Behav. Genet . 47 : 449 – 458 . OpenUrl PubMed ↵ Scheffers BR et al. 2016 . The broad footprint of climate change from genes to biomes to people . Science . 354 : aaf7671 . OpenUrl Abstract / FREE Full Text ↵ Scott JG . 1999 . Cytochromes P450 and insecticide resistance . Insect Biochem. Mol. Biol . 29 : 757 – 777 . OpenUrl CrossRef PubMed Web of Science ↵ Shahriari M , Zibaee A , Sahebzadeh N , Shamakhi L . 2018 . Effects of α-pinene, trans-anethole, and thymol as the essential oil constituents on antioxidant system and acetylcholine esterase of Ephestia kuehniella Zeller (Lepidoptera: Pyralidae) . Pestic. Biochem. Physiol . 150 : 40 – 47 . OpenUrl PubMed ↵ Sømme L . 1989 . Adaptations of terrestrial arthropods to the alpine environment . Biological Reviews . 64 : 367 – 407 . OpenUrl CrossRef ↵ Tan A , Abecasis GR , Kang HM . 2015 . Unified representation of genetic variants . Bioinformatics . 31 : 2202 – 2204 . OpenUrl CrossRef PubMed ↵ Tarasov A , Vilella AJ , Cuppen E , Nijman IJ , Prins P . 2015 . Sambamba: fast processing of NGS alignment formats . Bioinformatics . 31 : 2032 – 2034 . OpenUrl CrossRef PubMed ↵ Thomas CD et al. 2004 . Extinction risk from climate change . Nature . 427 : 145 – 148 . OpenUrl CrossRef PubMed Web of Science ↵ Torres-Banda V et al. 2022 . Gut transcriptome of two bark beetle species stimulated with the same kairomones reveals molecular differences in detoxification pathways . Comput. Struct. Biotechnol. J . 20 : 3080 – 3095 . OpenUrl PubMed ↵ Tuttle AM , Hoffman TL , Schilling TF . 2014 . Rabconnectin-3a regulates vesicle endocytosis and canonical Wnt signaling in zebrafish neural crest migration . PLoS Biol . 12 : e1001852 . OpenUrl CrossRef PubMed ↵ Vasseur P et al. 2022 . Human exposure to larvae of processionary moths in France: study of symptomatic cases registered by the French poison control centres between 2012 and 2019 . Clin. Toxicol . 60 : 231 – 238 . OpenUrl View the discussion thread. Back to top Previous Next Posted November 25, 2025. Download PDF Supplementary Material Email Thank you for your interest in spreading the word about bioRxiv. NOTE: Your email address is requested solely to identify you as the sender of this article. Your Email * Your Name * Send To * Enter multiple addresses on separate lines or separate them with commas. You are going to email the following A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth Message Subject (Your Name) has forwarded a page to you from bioRxiv Message Body (Your Name) thought you would like to see this page from the bioRxiv website. Your Personal Message CAPTCHA This question is for testing whether or not you are a human visitor and to prevent automated spam submissions. Share A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth Pierre Nouhaud , Mathieu Gautier , Alexandre Schifano , Laure Sauné , Carole Kerdelhué , Charles perrier bioRxiv 2025.11.24.690095; doi: https://doi.org/10.1101/2025.11.24.690095 Share This Article: Copy Citation Tools A major selective sweep likely linked to insecticide resistance identified along altitudinal gradients in the pine processionary moth Pierre Nouhaud , Mathieu Gautier , Alexandre Schifano , Laure Sauné , Carole Kerdelhué , Charles perrier bioRxiv 2025.11.24.690095; doi: https://doi.org/10.1101/2025.11.24.690095 Citation Manager Formats BibTeX Bookends EasyBib EndNote (tagged) EndNote 8 (xml) Medlars Mendeley Papers RefWorks Tagged Ref Manager RIS Zotero Tweet Widget Facebook Like Google Plus One Subject Area Evolutionary Biology Subject Areas All Articles Animal Behavior and Cognition (7635) Biochemistry (17690) Bioengineering (13892) Bioinformatics (41936) Biophysics (21451) Cancer Biology (18588) Cell Biology (25499) Clinical Trials (138) Developmental Biology (13378) Ecology (19899) Epidemiology (2067) Evolutionary Biology (24320) Genetics (15609) Genomics (22506) Immunology (17736) Microbiology (40394) Molecular Biology (17181) Neuroscience (88603) Paleontology (666) Pathology (2832) Pharmacology and Toxicology (4824) Physiology (7641) Plant Biology (15152) Scientific Communication and Education (2045) Synthetic Biology (4294) Systems Biology (9825) Zoology (2271)