Ecological Niche Adaptations Influence Transposable Element Dynamics in Pollinating and Non‐Pollinating Fig Wasps

Abstract

ABSTRACT

This study explores how ecological niches influence the dynamics of transposable elements (TEs) in the genomes of pollinating and non‐pollinating fig wasps (NPFWs), and how these ecological factors shape genome evolution. To examine the protective role of fig fruits for pollinators, we compared TE load and dynamics in six pollinating and five NPFW species from six different Ficus species. Phylogenetic analysis was used to assess correlations between genome size, oviposition sites, and TE length. We also analyzed the effects of natural selection and population dynamics on TE accumulation. Significant differences were observed in the total length, number, and types of TEs between pollinators and NPFWs. Phylogenetic analysis indicates that TEs in NPFWs, driven by genome size and oviposition sites, exhibit an expanding state, while pollinators show “dormant” TE landscapes with limited insertions. Despite relaxed selection pressure aimed at prolonging TE retention, pollinators maintain a limited TE abundance, likely due to the contracted population size. Additionally, numerous cis‐regulatory modules derived from TEs are located near genes involved in environmental information processing, emphasizing their potential role in adaptation. Our findings highlight the role of ecological niches, represented by oviposition sites, in shaping the TE dynamics of fig wasps. These results provide new insights into how ecological pressures influence genome evolution and adaptation in insects.

Introduction

Ecological niches play a critical role in shaping species' adaptation and evolutionary processes. However, how ecological niches influence genome evolution, particularly through transposable elements (TEs) and their temporal and spatial organization in the genome landscape, remains underexplored. TEs contribute significantly to genetic diversity and adaptation, and their activity and genome distribution can reveal not only how ecological pressures influence TE content but also the architecture of genomes. For example, in the invasive ant Cardiocondyla obscurior, bursts of TEs (TE islands) associated with olfactory receptor genes facilitate adaptation to new environments (Schrader et al. 2014). Similarly, in the aquatic plant Nymphoides indica, more than half of the genome consists of TEs, which may reflect adaptive responses to changing environmental conditions (Yang et al. 2022). In contrast, extreme environments, such as Antarctic ecosystems and enclosed fig fruits, are associated with compressed TE landscapes characterized by low diversity and content, as observed in the midge Belgica antarctica (Kelley et al. 2014) and the obligate mutualistic fig wasp Ceratosolen solmsi (Xiao et al. 2013). Notably, the specialized ecology of Drosophila incompta does not appear to significantly alter either the proportion and diversity of TEs in its genome compared to generalist Drosophila species, nor does it establish a distinct TE landscape pattern (Fonseca et al. 2019). These studies suggest that ecological pressures shape TE dynamics in different ways, although results across species have been inconsistent.

In insects, TE content and activity vary widely, leading to dynamic TE landscapes across species (Gilbert et al. 2021). These variations are shaped by a complex balance between TE acquisition, amplification, and loss influenced by both vertical inheritance and horizontal transfer (Walsh et al. 2013; Venner et al. 2017). Extensive studies in arthropods and insects indicate a robust phylogenetic signal and a positive correlation between genome size and TE content, supporting the concept of TEs as hereditary components transferred from parent to offspring (Sessegolo et al. 2016; Petersen et al. 2019; Wu and Lu 2019). However, a study on various Drosophila species suggests that at least one-third of TE families were acquired through horizontal transfer (Bartolome et al. 2009). The replication dynamics of TEs in the genome are intimately linked to their activity, with host-related factors such as immunity, natural selection, and genetic drift influencing the activity and fate (e.g., loss) of sequences inserted through TE amplification (Szitenberg et al. 2016; Arkhipova 2018; Kofler 2019; Jiang et al. 2024). The interplay of these processes reflects not only intrinsic genetic factors but also the host's ecological and evolutionary context.

The fig-fig wasp mutualism provides a unique system to explore the relationship between ecological niche and TE dynamics. Fig wasps are divided into two main groups: pollinating fig wasps (pollinators) and non-pollinating fig wasps (NPFWs) (Figure 1). As obligate mutualists, pollinators complete their entire life cycle within the syconium and serve as primary colonizers during the female flowering phase. Their morphological adaptations—flattened heads and robust hind legs—facilitate entry through the ostiole, while specialized pollen combs and baskets ensure effective pollination. In contrast, NPFWs exhibit more complex life cycles with broader host interactions. Most NPFWs cannot enter the syconium and have instead evolved elongated ovipositors with interspecific length variation to pierce the fig wall. This ecological divergence extends to niche differentiation: while pollinators colonize syconia synchronously with floral receptivity, NPFWs display distinct temporal and spatial oviposition strategies across ecological guilds (Cook and Rasplus 2003; Ranganathan et al. 2010; Wang et al. 2010). Previous studies have shown distinct genomic adaptations in these groups, with pollinators exhibiting reduced gene families for environmental sensing and detoxification, likely due to their stable and protected habitat inside figs (Xiao et al. 2021). Given these contrasting ecological niches, we hypothesize that: (i) TE landscape divergence. Pollinators, living in a closed, stable environment, will have a “dormant” TE landscape, characterized by the predominance of ancient TE remnants and limited recent TE expansion. In contrast, NPFWs with external oviposition strategies that encounter more variable external environments are expected to exhibit an “aggressive” TE landscape, with ongoing TE insertions and stress-responsive TE activation. (ii) Circadian system specialization. We further postulate that pollinators' strict oviposition timing—requiring precise synchronization with the brief female phase (typically < 48 h)—has been selected for enhanced circadian clock gene conservation and reduced flexibility in activity rhythms compared to generalist NPFWs.

[IMAGE OMITTED. SEE PDF]

This study aims to analyze the composition, landscape characteristics, and influencing factors of TEs in the genomes of pollinating and NPFWs. By comparing the TE accumulation and activity between these two groups, we seek to clarify the role of spatial oviposition niches (oviposition sites) in shaping TE dynamics. We also explore how pollinators' precise phenological synchronization relates to the conservation of circadian clock genes and TE regulation mechanisms. Our results will provide new insights into how ecological pressures influence genome evolution and species adaptation, laying the groundwork for future studies on biogeography and conservation in changing environments.

Materials and Methods

Sample Collection and Genome Assembly

We analyzed six species of pollinating fig wasps (pollinators) and five species of NPFWs, representing all major feeding habits, diverse host fig species, and contrasting oviposition sites. Details of the selected species are presented in Table 1. The voucher specimens have been deposited in the Institute of Entomology, College of Life Sciences, Nankai University, Tianjin, China.

TABLE 1 Details of fig wasps and TEs annotations.

Category	Species	Abbreviation	Host fig tree	Collection locality	Collection date	Mating site	Oviposition site	Feeding habit	r	μ	Genome size (Mb)	TE length (Mb)	No. of TEs	No. of classified TE families	TE content (%)
Pollinators	Platyscapa corneri	P. corneri	Ficus superba	Mount Bawang, Hainan	2016/7/11	Inside	Inside	Herbivorous	0.0049	2.43E-09	304.37	18.23	57,619	56	5.99
Eupristina koningsbergeri	E. koningsbergeri	Ficus benjamina	Danzhou, Hainan	2016/9/14	Inside	Inside	Herbivorous	0.0051	2.55E-09	318.36	48.99	113,791	54	15.39
Wiebesia pumilae	W. pumilae	Ficus pumila	Huangshan, Anhui	2017/4/10	Inside	Inside	Herbivorous	0.0047	2.37E-09	319.91	25.07	80,474	55	7.84
Dolichoris vasculosae	D. vasculosae	Ficus vasculosa	Mount Jianfeng, Hainan	2017/4/3	Inside	Inside	Herbivorous	0.0050	2.48E-09	286.21	8.02	30,895	50	2.80
Kradibia gibbosae	K. gibbosae	Ficus gibbosa	Xishuangbanna, Yunnan	2016/12/1	Inside	Inside	Herbivorous	0.0052	2.60E-09	230.30	7.66	42,828	53	3.33
Ceratosolen fusciceps	C. fusciceps	Ficus racemosa	Xishuangbanna, Yunnan	2016/12/19	Inside	Inside	Herbivorous	0.0053	2.66E-09	235.18	9.06	57,250	52	3.85
NPFWs	Sycobia sp.2	Sycobia sp.2	Ficus benjamina	Nada, Hainan	2017/5/18	Inside/Outside	Outside	Herbivorous	0.0048	2.40E-09	621.39	235.27	673,913	63	37.86
Philotrypesis tridentata	P. tridentata	Ficus benjamina	Danzhou, Hainan	2017/5/4	Inside/Outside	Outside	Omnivorous	0.0049	2.45E-09	398.61	170.35	401,464	65	42.74
Apocrypta bakeri	A. bakeri	Ficus hispida	Xinglong, Hainan	2017/5/6	Inside	Outside	Omnivorous	0.0051	2.57E-09	198.59	27.51	146,704	61	13.85
Sycophaga agraensis	S. agraensis	Ficus racemosa	Xishuangbanna, Yunnan	2017/7/1	Inside	Outside	Carnivorous	0.0048	2.40E-09	246.14	45.80	114,334	60	18.61
Sycophila sp.2	Sycophila sp.2	Ficus benjamina	Nada, Hainan	2017/6/10	Inside/Outside	Outside	Omnivorous	0.0043	2.16E-09	288.59	80.51	252,000	63	27.90

The whole bodies of approximately 150 female fig wasps were pooled for DNA extraction for each species. The short-insert paired-end library (450 bp) was constructed using the Illumina Kit (San Diego, CA, USA) and sequenced on an Illumina HiSeq 2500 platform, generating ~50× coverage. The long-insert library (20 kb) was prepared using SMRTbell Express Template Prep Kits (Pacific Biosciences, Menlo Park, CA, USA) and run on a PacBio Sequel system, yielding ~50× coverage.

The PacBio data were corrected with Illumina data by proovread (v2.14.0) (Hackl et al. 2014) and then assembled using SMRTdenovo (v1.0) (). Two polishing strategies were performed with Arrow (Chin et al. 2013) and Pilon (v1.22) (Walker et al. 2014) to correct the assembly with PacBio long reads and Illumina short reads, respectively. The quality and completeness of each genome assembly were assessed by comparing it with BUSCO (v2.0) (Manni et al. 2021) genes, achieving over 95.9% completeness with more than 94.0% single-copy genes. Duplicated BUSCO ratios ranged from 0.8% to 1.9%.

Transposable Element Annotation

TEs were annotated by constructing de novo repeat libraries using RepeatModeler (, v2.0.2) with the “-LTRStruct” parameter. Libraries were filtered using species-specific protein annotation files and merged with the “chalcidoidea” TE library, extracted using the Perl script queryRepeatDatabase.pl. from RepeatMasker's util file. We then identified TE families by combining these libraries with RepBase (v20181026) and Dfam (v3.5) databases. For accuracy, the genomes were masked with RepeatMasker (v4.1.2-p1), using a cutoff parameter of 250. TE annotation results were summarized with the parseRM.pl. () script. Differences in total length, number, types, and genomic content of TEs between pollinators and NPFWs were statistically evaluated using the Mann–Whitney U test.

TE Insertion Time Estimation

The insertion times (T) were estimated using the formula T = K/2r (Li et al. 2019), where divergence (K) was calculated from Kimura's formula K = –300/4 × Ln(1 – D × 4/300), with mutation rates (perc div.) derived for each TE sequence. The nucleotide substitution rate (r) was calculated based on 4-fold degenerate sites from single-copy orthologs of 17 species, as described by Xiao et al. (2021). Branch lengths were estimated using PhyloFit in PHAST (, v1.5), and the root-to-tip branch length was determined with TreeStat (, v1.2). The nucleotide substitution rate was then computed by dividing the root-to-tip branch length by the divergence time of D. pulex (485 Mya) on the phylogenetic tree (Xiao et al. 2021).

Phylogenetic Analyses

The phylogenetic relationship of 12 fig wasps (including 11 fig wasps in this research and Ceratosolen solmsi) was reconstructed using single-copy orthologous genes, with Apis mellifera designated as the outgroup. Single-copy orthologs were identified using OrthoMCL (v2.0.9) (Li et al. 2003). Protein sequences of each gene were aligned using MAFFT (v7.313) (Katoh and Standley 2013) and poorly aligned regions were removed using Gblocks (v0.91b) (Talavera and Castresana 2007). A maximum likelihood tree was constructed using IQ-TREE (v1.6.12) (Nguyen et al. 2015) with JTT + F + I + G4 model and 5000 bootstrap replicates.

A phylogenetic tree used for phylogenetic comparative analysis was constructed using TreeViewer (v2.2.0) (Bianchini and Sanchez-Baracaldo 2024). The phylogenetic signal of TE length, genome size, and oviposition sites was evaluated using Blomberg's K and D values (Fritz and Purvis 2010). A phylogenetic tree-based Monte Carlo simulation approach was used to evaluate the statistical power of Blomberg's K estimation. Using the fastBM function from the phytools (Revell 2024) package, we generated continuous trait data under Brownian motion models with specified true K values (0.9, 1.0, and 1.1). For each simulated dataset, phylogenetic signal strength was assessed via phylosig with significance threshold alpha = 0.05 (testing H₀: K = 0). After 1000 simulation replicates, statistical power was quantified as the proportion of iterations rejecting H₀.

Relationships between TE length, genome size, and oviposition sites were analyzed using phylogenetic generalized least squares (PGLS) in R with the “APE” (Paradis and Schliep 2019), “phytools” (Revell 2024), and “caper” (Freckleton et al. 2002) packages. Model selection was based on the Akaike Information Criterion corrected (AICc) for small sample sizes. Ancestral trait reconstruction of TE length was performed using the results predicted by the multifactor PGLS model.

Natural Selection Analysis

To detect selection differences between pollinators and NPFWs, we used Codeml from PAML (v4.10.7) (Yang 2007) and the RELAX model in HyPhy (v2.5.57) (Pond et al. 2020). Single-copy orthologous genes from the 11 species were aligned using MAFFT (v7.520) (Katoh and Standley 2013) and filtered with Gblocks (v0.91b) (Talavera and Castresana 2007) to remove poorly aligned regions. Codeml's free-ratio model was used to estimate nonsynonymous-to-synonymous substitution ratios (dN/dS ratios), and differences between pollinators and NPFWs were compared using the Mann–Whitney U test.

Positive selection genes in pollinator lineages were identified using the branch model in Codeml through the likelihood ratio test. Selection intensity analysis was conducted on pollinator branches with the RELAX model, and gene sets under relaxed and intensified selection were collected. To compare the effects of selection pressure on TEs, the lengths of TEs within 1 Kb of these genes were analyzed. To account for differences in TE types near different genes, TE lengths were normalized using consensus sequences.

Demographic History Analysis

Demographic histories were reconstructed for each species using pairwise sequential markovian coalescence (PSMC, v0.6.5) (Li and Durbin 2011) based on whole-genome resequencing data. Genomic DNA was extracted using the DNeasy blood and tissue kit (QIAGEN GmbH, QIAGEN Strasse 1, 40 724 Hilden, GERMANY) and sequenced with DNBSEQ (DNBSEQ platform developed by BGI) by BGI Genomics Co. Ltd. Reads were aligned to the genome using BWA-MEM (v0.7.17) (Li and Durbin 2009), and SNP variants were called with BCFtools (v1.8) (Danecek et al. 2021). PSMC was run with specific atomic intervals “4 + 25 × 2 + 4 + 6,” using the neutral mutation rate (μ) calculated from the nucleotide substitution rate per million years per site (r).

cis-Regulatory Module Prediction

cis-regulatory modules (CRMs) were predicted using SCRMshaw_HD (, vHD) (Kazemian and Halfon 2019) with tandem repeats masked by Tandem Repeats Finder (TRF, v4.09) (Benson 1999). Predictions were based on Drosophila training data, with window size (--wlen) set to 200 bp with 100 bp shifts (--wshift). To avoid overprediction, the maximum number of hits for top-scored windows (--thitw) was set to 2000. Finally, the needed hits chosen by Generate_top_N_SCRMhits.pl script were merged with BEDTools (v2.19.1) (Quinlan and Hall 2010) and compared to TE locations to assess overlapping patterns.

Gene Function Enrichment Analysis

For functional enrichment analysis, genes within 1 Kb of classified TEs from major burst peaks and CRMs overlapping with TEs were selected. GO and KEGG enrichment analyses were performed using ClueGO (Bindea et al. 2009) in Cytoscape (v3.9.1) with default parameters. Term-term interrelation was defined based on shared genes using the kappa score. For better visualization of commonalities, the enrichment networks were merged in the intersection mode, and the nodes were filtered using the parameter Term p-value < 0.05.

Results

TE Composition of Fig Wasps Across Different Ecological Niches

The total TE length in six pollinating fig wasp (pollinators) species varied between 7.66 and 48.99 Mb, while in five NPFWs species ranged from 27.51 to 235.27 Mb (Table 1). Mann–Whitney U test revealed that both the total length (U = 28.000, Z = 2.373, p = 0.017, Figure S1a) and the number (U = 30.000, Z = 2.739, p = 0.004, Figure S1b) of TEs were significantly higher in NPFWs than in pollinators. Pollinators exhibited between 50 and 56 classified TE families, whereas NPFWs ranged from 60 to 65, with significant differences between the two groups (Mann–Whitney U = 30.000, Z = 2.745, p = 0.004, Figure S1c). Among the classified TE families, 46 types were common across all species, while six types were exclusive to NPFWs. The Crypton-I family was absent in pollinators, except for one detected instance in Wiebesia pumilae after parameter adjustments, indicating a trend of TE loss in pollinators.

The genomic content of TEs in pollinators ranged from 2.80% to 15.39%, while in NPFWs, it was significantly higher, ranging from 13.85% to 42.74% (Table 1). The Mann–Whitney U test confirmed that these differences in TE content between the two groups were statistically significant (U = 29.000, Z = 2.556, p = 0.009, Figure S1d).

Within the Ficus benjamina system, the TE content of three NPFWs (Sycobia sp.2, P. tridentata, and Sycophila sp.2) (mean ± SD: 36.17% ± 7.56%) significantly exceeded that of their obligate pollinator (E. koningsbergeri: 15.39%) (one-sample t-test: t = 4.758, df = 2, p = 0.041; Table S1), with a 135% elevation in TE load (95% CI: 17.38%–54.95%). This host-controlled comparison suggests niche-specific (pollinator vs. NPFWs) TE accumulation independent of host phylogeny.

Ancestral TE Load Reconstruction Indicated Slower Increase in Pollinators

Phylogenetic signal analysis revealed distinct evolutionary patterns across traits. Neither TE length nor genome size showed a significant phylogenetic signal (Blomberg's K for TE length = 0.948, p = 0.176; for genome size = 0.863, p = 0.357). Phylogenetic tree-based Monte Carlo simulations revealed limited statistical power (27.4%–29.2%) for detecting Blomberg's K values in the range of 0.9–1.1. However, a significant phylogenetic signal was observed for oviposition sites (D = −5.526, P_random = 0, P_Brownian = 0.977).

Using various PGLS models to explore their relationships, we found that genome size and oviposition sites were independently correlated with TE length (genome size: λ = 1.00, F = 137.10, p < 0.001; oviposition sites: λ = 0.00, F = 26.06, p < 0.01), whereas their interaction was not significant (λ = 0.00; F = 1.23, p = 0.303). The lowest AICc value (97.34) of the multifactor model confirmed the independent effects of genome size and oviposition sites on TE length. The coefficient for genome size (0.495) was notably smaller than that for the “outside” oviposition site (58.594), suggesting that the categorical effect of the “outside” mode had a stronger influence on TE length. As shown in the GLS model (Figure 2C), the increase in TE length in the NPFWs group was faster with the increase in genome size compared to pollinators. Ancestral state reconstruction of TE length (Figure 2A,D), based on the PGLS multifactor model, indicated that TE length in pollinators was shorter than in NPFWs.

[IMAGE OMITTED. SEE PDF]

Active TEs Shape Distinct Landscapes

Distinct temporal patterns emerged between pollinators and NPFWs (Figure 3, outside panel). Pollinators exhibited ancient “dormant” peaks (median 77.91% of total TE length from 10 to 30 Mya), reflecting early TEs invasion and amplification followed by suppression and decay, with limited recent TE insertions. In contrast, NPFWs displayed recent “aggressive” peaks, with 17.68% of TE load originating from the last 5 million years (vs. 8.18% in pollinators), highlighting a periodically recurrent influx of active TE insertions.

[IMAGE OMITTED. SEE PDF]

The dominant TE types varied among species, yet a consistent pattern emerged where 2–4 TE families collectively accounted for over 50% of classified TEs in all species (Table S2). For example, P. corneri exhibited a high abundance of the RTE-X family, with a total length of 4.97 Mb, accounting for 42.7% of classified TEs, while D. vasculosae and K. gibbosae showed dominance of the I-Jockey and R1-LOA families, respectively. Notably, in all NPFWs and several pollinators, the Gypsy family consistently accounted for the highest content, suggesting that bursts of specific TE types have shaped each species' TE landscape (Figure 3, inside panel). Furthermore, landscapes of Gypsy family members in NPFWs (Figure S2) showed that a few members have undergone recent bursts (Figure S2 outside panels), while the multiple bursts of individual members over time indicated sustained genomic activity (Figure S2 inside panels).

Influence of Natural Selection on TEs Accumulation in Pollinators

Negative selection was the dominant evolutionary force shaping TE-related genes in both pollinators and NPFWs (Figure S3). However, nonsynonymous-to-synonymous substitution ratios (dN/dS) were significantly higher in pollinators than in NPFWs (ANOVA, F = 30.365, p < 0.001) (Figure S3). Subsequently, gene sets under relaxed selection (1129), intensified selection (512), and positive selection (642) were identified from 4238 single-copy homologous genes. The lengths of TEs within 1Kb proximity to these gene sets were compared. Multiple comparisons following ANOVA (Figure 4) revealed that except for D. vasculosae, TEs near genes under relaxed selection exhibited significantly longer lengths compared to those near genes under intensified or positive selection. Generally, TEs near positively selected genes had the shortest length, with the exception of K. gibbosae.

[IMAGE OMITTED. SEE PDF]

Historical Population Dynamics of Fig Wasps

Effective population size (Ne) trajectories revealed distinct evolutionary histories between pollinators and NPFWs over the past 1 million years (Mya) (Figure 5). Pollinators have generally exhibited a decline in Ne, with two distinct patterns: a direct decline, as seen in P. corneri, and a fluctuating decline, as in W. pumilae. In contrast, the Ne trajectories of NPFWs were more complex. Some species, such as P. tridentata, A. bakeri, and S. agraensis, showed increases amidst fluctuations, while others, like Sycobia sp.2 and Sycophila sp.2, displayed fluctuating declines. Recent relative Ne values for pollinators ranged from 1.306 × 10⁴ to 9.626 × 10⁴, significantly lower than those for NPFWs, which ranged from 8.717 × 10⁴ to 46.919 × 10⁴ (Mann–Whitney U = 29.000, Z = 2.556, p = 0.009).

[IMAGE OMITTED. SEE PDF]

Functional Analysis of Genes Near Burst TEs in Pollinators

Temporal patterns of TE burst varied among pollinators, with peak activities at 24–26 Mya (P. corneri), 20–21 Mya (E. koningsbergeri), 20–22 Mya (W. pumilae), 11–12 Mya (D. vasculosae), 20–21 Mya (K. gibbosae), and 15–16 Mya (C. fusciceps). Genes adjacent to TE insertion hotspots exhibited significant functional enrichment across pollinator species. Gene Ontology (GO) enrichment analysis revealed significant enrichment in pathways related to ion transport, signal transduction, etc. Similarly, the Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis indicated significant enrichment in pathways associated with environmental information processing, such as circadian entrainment (term ID: 04713) (Figure S4).

Notably, among the genes enriched in the circadian entrainment pathway, the NOS1 (Nitric Oxide Synthase 1) gene was shared by 5 out of 6 pollinator species. Domain architecture analysis of NOS1 (Table S4) revealed striking differences between pollinators and NPFWs. In pollinators, the NOS1 domains were highly conserved: Domain 1 consistently belonged to the NOS_oxygenase superfamily (cl00254 or cd00795), Domain 2 (FNR_like, flavin-binding) was retained in all species except for E. koningsbergeri, and Domain 3 (electron transfer) featured either CysJ (cl43121) or Flavodoxin_1 (pfam00258). In contrast, NPFWs exhibited frequent domain loss or divergence.

The Evolutionary Relationship Between TEs and Cis-Regulatory Modules

Substantial overlap was observed between TEs and CRMs, with 5.4%–14.1% of TEs intersecting CRMs across pollinator species (Table S3). Kradibia gibbosae exhibited the highest proportion of overlap, and its genome-wide distribution revealed a reciprocal density pattern: regions with high CRM density tended to have fewer TEs, and vice versa. Despite this pattern, substantial overlap was observed in three configurations (Figure 6). The most common configuration was TEs located within CRMs, accounting for nearly 60% of cases, followed by regions where TEs and CRMs intersected. The least prevalent configuration was CRMs within TEs, representing only 11.9%. Furthermore, gene function enrichment analysis of regions near CRMs overlapping with TEs consistently highlighted pathways involved in environmental information processing.

[IMAGE OMITTED. SEE PDF]

Discussion

Our analysis of TEs in fig wasps across different ecological niches provides new insights into the factors shaping TE dynamics, including genome size, oviposition sites, selective pressures, and effective population size. The observed differences in TE load and landscape between pollinators and NPFWs highlight the profound influence of ecological constraints and evolutionary histories on genome architecture.

Ecological Niche and TE Dynamics

The lower TE content and family diversity in pollinators, along with the absence of Crypton-I, may reflect genomic streamlining due to strict vertical transmission and co-evolution with their host figs. Although Mann–Whitney U tests have limitations when analyzing phylogenetically nonindependent, non-monophyletic NPFW species, our host-controlled comparison within the F. benjamina system effectively minimized phylogenetic bias while strongly supporting the “niche effect” hypothesis. Despite current sample size constraints, this case study offers compelling preliminary evidence that TE divergence can arise rapidly within shared host environments. These findings are consistent with patterns observed in other host-specialized insect systems, where niche-specific selection pressures appear to supersede phylogenetic history (Jousselin et al. 2013).

Our phylogenetic analysis revealed intriguing evolutionary patterns in fig wasp traits. While neither TE length nor genome size exhibited a significant phylogenetic signal, this pattern requires careful interpretation. Monte Carlo simulations confirmed our analysis had limited power to detect moderate phylogenetic signals, suggesting these results may reflect either genuine evolutionary lability or methodological constraints inherent to small clade studies. Despite this limitation, our PGLS models suggest that TE length is still influenced by genome size in a way that aligns with phylogeny. This finding supports the notion that larger genomes tend to accommodate more TEs (Elliott and Gregory 2015). However, the most significant predictor of TE length was the oviposition sites, with external oviposition (“outside” mode) linked to larger TEs. This suggests that exposure to external environmental pressures may facilitate or maintain higher TE activity (Roquis et al. 2021), supporting the hypothesis that TE dynamics are shaped by the complexity and stability of ecological niches.

Pollinators, living in a stable and closed environment within fig fruits, may be subject to stronger suppression mechanisms due to the limited need for rapid adaptation (Bourque et al. 2018), and exhibit a “dormant” TE landscape, characterized by ancient TE insertions that have likely decayed over time. In contrast, NPFWs face diverse parasitic strategies and complex feeding habits, with some oligophagous species interacting with more than one host fig species, necessitating heightened abilities for host and food source detection (Borges 2015). Moreover, the extended duration that NPFWs spend on the external surfaces of syconia for mating and oviposition exposes them to predation and disease risks, potentially driving a need for enhanced defense mechanisms (Ranganathan et al. 2010). Additionally, gallers such as pollinators are proovigenic—emerging with their full complement of mature eggs—while parasitoids are synovigenic, maturing eggs gradully over time (Borges 2015). This difference may be related to the extent to which their eggs are exposed to external environment influences. In addition to the factors discussed above, other ecological traits—particularly the extended adult lifespan and foraging requirements of NPFWs—may further contribute to observed TE dynamics. It has been reported that TEs become active under stress conditions (Horvath et al. 2017), and the active state is a prerequisite for the burst of TEs (Belyayev 2014). This is consistent with the “aggressive” TE landscape, with recent bursts and ongoing TE activity. For instance, the Gypsy family of TEs, which accounts for a large proportion of the total TE load in NPFWs, has undergone multiple bursts, highlighting the active nature of these elements. This indicates that TEs in NPFWs are more likely to be driven by intrinsic TE activity, rather than horizontal transfer.

Selection Pressures and TE Accumulation

Selection pressures, especially purifying selection, often limit TE accumulation by removing deleterious or neutral TEs from the genome (Stritt et al. 2018). Our analysis shows that both groups are predominantly under negative selection, with pollinators exhibiting higher dN/dS ratios, which means more relaxed purifying selection and indeed prolonged retention of TEs within their genomes. This trend indicates that TE accumulation is modulated by the selective context of neighboring genes (Baduel et al. 2019). However, it was apparent that this factor alone did not serve as the explanation for the observed disparity in TE load between the two distinct categories of fig wasps.

Role of Effective Population Size

Effective population size also has significant and multifaceted effects on the accumulation of genetic variation. In smaller populations, the strength of genetic drift is higher, which can obscure the effects of natural selection and allow deleterious mutations to persist (Mathur et al. 2023). However, the smaller Ne combined with relaxed purifying selection clearly cannot explain the lower TE load in pollinators. An alternative perspective suggests that when populations contract, deleterious mutations are more likely to be eliminated (Xie et al. 2022). This is more likely the reason for the lower TE load in pollinators.

TEs as Potential Drivers of Adaptation

The concentrated peak of residual TEs in the genomes of pollinator groups over a specific timescale piqued our curiosity. Considering that TE insertions or deletions could emerge as a potent factor impacting the expression of neighboring genes (Hollister and Gaut 2009; Oliver and Greene 2011), we conducted an in-depth analysis of TE fragments concentrated from approximately 12 to 24 million years ago (Mya) in pollinators. A predominant association was observed between these TEs and genes involved in signal transduction and environmental information processing. Further investigation into the relationship between TEs and CRMs suggests that these residual TEs may function as CRMs, influencing the expression of nearby genes in response to environmental stress. An illustrative example is the role of TEs in circadian entrainment; as climate-induced changes affect fig phenology, circadian rhythm regulation may help fig wasps adapt to these shifts (Wang et al. 2024). The highly conserved NOS1 gene in pollinators serves as compelling evidence when compared to NPFWs. This supports the idea that TEs, far from being mere genomic “junk,” may serve important regulatory functions that contribute to host adaptation in response to ecological niche pressures.

Conclusion

This study underscores the significant impact of ecological niches on the accumulation of TEs in fig wasps. The correlation between genome size and TE content reflects a universal genomic scaling law, but in fig wasps, this relationship seemingly functions independently of phylogenetic constraints, an aspect that warrants further investigation. Instead, the variations in TE activity are primarily driven by ecological pressures and niche-specific factors. These findings reveal how ecological differences influence the complex dynamics of TE accumulation. Future research could refine methods to track TE gains and losses over time, providing deeper insights into TE evolution and its role in host adaptation.

Author Contributions

Jing Liu: data curation (equal), formal analysis (lead), funding acquisition (lead), investigation (equal), methodology (lead), software (equal), visualization (lead), writing – original draft (lead), writing – review and editing (lead). Yun-Heng Miao: formal analysis (supporting), methodology (supporting), software (equal). Hong-Xia Hou: formal analysis (supporting), funding acquisition (supporting), investigation (equal). Da-Wei Huang: funding acquisition (equal), project administration (lead), writing – review and editing (equal). Jin-Hua Xiao: project administration (lead), resources (equal), writing – review and editing (lead).

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant numbers: 32301412, 32070466, 32200375) and the Hebei Natural Science Foundation (Grant number: C2024108006).

Conflicts of Interest

The authors declare no conflicts of interest.

Data Availability Statement

Data are available under NCBI BioProjects PRJNA641212 and PRJNA494992. Currently released files include historical genome assemblies; manuscript-specific data will be fully accessible upon publication (or to 2027-12-31). Accession numbers: SAMN44830401, SAMN44830400, SAMN44830399, SAMN44830398, SAMN44830396, SAMN44830397, SAMN44830395, SAMN44830402, SAMN44830394, SAMN44830393, SAMN44830392.

References

Arkhipova, I. R. 2018. “Neutral Theory, Transposable Elements, and Eukaryotic Genome Evolution.” Molecular Biology and Evolution 35, no. 6: 1332–1337. https://doi.org/10.1093/molbev/msy083.

Baduel, P., L. Quadrana, B. Hunter, K. Bomblies, and V. Colot. 2019. “Relaxed Purifying Selection in Autopolyploids Drives Transposable Element Over‐Accumulation Which Provides Variants for Local Adaptation.” Nature Communications 10, no. 1: 5818. https://doi.org/10.1038/s41467‐019‐13730‐0.

Bartolome, C., X. Bello, and X. Maside. 2009. “Widespread Evidence for Horizontal Transfer of Transposable Elements Across Drosophila Genomes.” Genome Biology 10, no. 2: R22. https://doi.org/10.1186/gb‐2009‐10‐2‐r22.

Belyayev, A. 2014. “Bursts of Transposable Elements as an Evolutionary Driving Force.” Journal of Evolutionary Biology 27, no. 12: 2573–2584. https://doi.org/10.1111/jeb.12513.

Benson, G. 1999. “Tandem Repeats Finder: A Program to Analyze DNA Sequences.” Nucleic Acids Research 27, no. 2: 573–580. https://doi.org/10.1093/nar/27.2.573.

Bianchini, G., and P. Sanchez‐Baracaldo. 2024. “TreeViewer: Flexible, Modular Software to Visualise and Manipulate Phylogenetic Trees.” Ecology and Evolution 14, no. 2: e10873. https://doi.org/10.1002/ece3.10873.

Bindea, G., B. Mlecnik, H. Hackl, et al. 2009. “ClueGO: A Cytoscape Plug‐In to Decipher Functionally Grouped Gene Ontology and Pathway Annotation Networks.” Bioinformatics 25, no. 8: 1091–1093. https://doi.org/10.1093/bioinformatics/btp101.

Borges, R. M. 2015. “How to Be a Fig Wasp Parasite on the Fig‐Fig Wasp Mutualism.” Current Opinion in Insect Science 8: 34–40. https://doi.org/10.1016/j.cois.2015.01.011.

Bourque, G., K. H. Burns, M. Gehring, et al. 2018. “Ten Things You Should Know About Transposable Elements.” Genome Biology 19, no. 1: 199. https://doi.org/10.1186/s13059‐018‐1577‐z.

Chin, C. S., D. H. Alexander, P. Marks, et al. 2013. “Nonhybrid, Finished Microbial Genome Assemblies From Long‐Read SMRT Sequencing Data.” Nature Methods 10, no. 6: 563. https://doi.org/10.1038/Nmeth.2474.

Cook, J. M., and J. Y. Rasplus. 2003. “Mutualists With Attitude: Coevolving Fig Wasps and Figs.” Trends in Ecology & Evolution 18, no. 5: 241–248. https://doi.org/10.1016/S0169‐5347(03)00062‐4.

Danecek, P., J. K. Bonfield, J. Liddle, et al. 2021. “Twelve Years of SAMtools and BCFtools.” GigaScience 10, no. 2: 1–4. https://doi.org/10.1093/gigascience/giab008.

Elliott, T. A., and T. R. Gregory. 2015. “Do Larger Genomes Contain More Diverse Transposable Elements?” BMC Evolutionary Biology 15: 69. https://doi.org/10.1186/s12862‐015‐0339‐8.

Fonseca, P. M., R. D. Moura, G. L. Wallau, and E. L. S. Loreto. 2019. “The Mobilome of Drosophila Incompta, a Flower‐Breeding Species: Comparison of Transposable Element Landscapes Among Generalist and Specialist Flies.” Chromosome Research 27, no. 3: 203–219. https://doi.org/10.1007/s10577‐019‐09609‐x.

Freckleton, R. P., P. H. Harvey, and M. Pagel. 2002. “Phylogenetic Analysis and Comparative Data: A Test and Review of Evidence.” American Naturalist 160, no. 6: 712–726. https://doi.org/10.1086/343873.

Fritz, S. A., and A. Purvis. 2010. “Phylogenetic Diversity Does Not Capture Body Size Variation at Risk in the World's Mammals.” Proceedings of the Biological Sciences 277, no. 1693: 2435–2441. https://doi.org/10.1098/rspb.2010.0030.

Gilbert, C., J. Peccoud, and R. Cordaux. 2021. “Transposable Elements and the Evolution of Insects.” Annual Review of Entomology 66: 355–372. https://doi.org/10.1146/annurev‐ento‐070720‐074650.

Hackl, T., R. Hedrich, J. Schultz, and F. Förster. 2014. “Proovread: Large‐Scale High‐Accuracy PacBio Correction Through Iterative Short Read Consensus.” Bioinformatics 30, no. 21: 3004–3011. https://doi.org/10.1093/bioinformatics/btu392.

Hollister, J. D., and B. S. Gaut. 2009. “Epigenetic Silencing of Transposable Elements: A Trade‐Off Between Reduced Transposition and Deleterious Effects on Neighboring Gene Expression.” Genome Research 19, no. 8: 1419–1428. https://doi.org/10.1101/gr.091678.109.

Horvath, V., M. Merenciano, and J. Gonzalez. 2017. “Revisiting the Relationship Between Transposable Elements and the Eukaryotic Stress Response.” Trends in Genetics 33, no. 11: 832–841. https://doi.org/10.1016/j.tig.2017.08.007.

Jiang, J., Y. C. Xu, Z. Q. Zhang, et al. 2024. “Forces Driving Transposable Element Load Variation During Arabidopsis Range Expansion.” Plant Cell 36, no. 4: 840–862. https://doi.org/10.1093/plcell/koad296.

Jousselin, E., A. Cruaud, G. Genson, F. Chevenet, R. G. Foottit, and A. Cœurd'acier. 2013. “Is Ecological Speciation a Major Trend in Aphids? Insights From a Molecular Phylogeny of the Conifer‐Feeding Genus.” Frontiers in Zoology 10: 56. https://doi.org/10.1186/1742‐9994‐10‐56.

Katoh, K., and D. M. Standley. 2013. “MAFFT Multiple Sequence Alignment Software Version 7: Improvements in Performance and Usability.” Molecular Biology and Evolution 30, no. 4: 772–780. https://doi.org/10.1093/molbev/mst010.

Kazemian, M., and M. S. Halfon. 2019. “CRM Discovery Beyond Model Insects.” Methods in Molecular Biology 1858: 117–139. https://doi.org/10.1007/978‐1‐4939‐8775‐7_10.

Kelley, J. L., J. T. Peyton, A. S. Fiston‐Lavier, et al. 2014. “Compact Genome of the Antarctic Midge Is Likely an Adaptation to an Extreme Environment.” Nature Communications 5: 4611. https://doi.org/10.1038/ncomms5611.

Kofler, R. 2019. “Dynamics of Transposable Element Invasions With piRNA Clusters.” Molecular Biology and Evolution 36, no. 7: 1457–1472. https://doi.org/10.1093/molbev/msz079.

Li, H., and R. Durbin. 2009. “Fast and Accurate Short Read Alignment With Burrows‐Wheeler Transform.” Bioinformatics 25, no. 14: 1754–1760. https://doi.org/10.1093/bioinformatics/btp324.

Li, H., and R. Durbin. 2011. “Inference of Human Population History From Individual Whole‐Genome Sequences.” Nature 475, no. 7357: 493–496. https://doi.org/10.1038/nature10231.

Li, L., C. J. Stoeckert Jr., and D. S. Roos. 2003. “OrthoMCL: Identification of Ortholog Groups for Eukaryotic Genomes.” Genome Research 13, no. 9: 2178–2189. https://doi.org/10.1101/gr.1224503.

Li, S. F., Y. J. Guo, J. R. Li, et al. 2019. “The Landscape of Transposable Elements and Satellite DNAs in the Genome of a Dioecious Plant Spinach (Spinacia Oleracea L.).” Mobile DNA 10: 3. https://doi.org/10.1186/s13100‐019‐0147‐6.

Manni, M., M. R. Berkeley, M. Seppey, F. A. Simao, and E. M. Zdobnov. 2021. “BUSCO Update: Novel and Streamlined Workflows Along With Broader and Deeper Phylogenetic Coverage for Scoring of Eukaryotic, Prokaryotic, and Viral Genomes.” Molecular Biology and Evolution 38, no. 10: 4647–4654. https://doi.org/10.1093/molbev/msab199.

Mathur, S., J. M. Tomecek, L. A. Tarango‐Arambula, R. M. Perez, and J. A. DeWoody. 2023. “An Evolutionary Perspective on Genetic Load in Small, Isolated Populations as Informed by Whole Genome Resequencing and Forward‐Time Simulations.” 77, no. 3: 690–704. https://doi.org/10.1093/evolut/qpac061.

Nguyen, L. T., H. A. Schmidt, A. von Haeseler, and B. Q. Minh. 2015. “IQ‐TREE: A Fast and Effective Stochastic Algorithm for Estimating Maximum‐Likelihood Phylogenies.” Molecular Biology and Evolution 32, no. 1: 268–274. https://doi.org/10.1093/molbev/msu300.

Oliver, K. R., and W. K. Greene. 2011. “Mobile DNA and the TE‐Thrust Hypothesis: Supporting Evidence From the Primates.” Mobile DNA 2, no. 1: 8. https://doi.org/10.1186/1759‐8753‐2‐8.

Paradis, E., and K. Schliep. 2019. “Ape 5.0: An Environment for Modern Phylogenetics and Evolutionary Analyses in R.” Bioinformatics 35, no. 3: 526–528. https://doi.org/10.1093/bioinformatics/bty633.

Petersen, M., D. Armisen, R. A. Gibbs, et al. 2019. “Diversity and Evolution of the Transposable Element Repertoire in Arthropods With Particular Reference to Insects.” BMC Evolutionary Biology 19, no. 1: 11. https://doi.org/10.1186/s12862‐018‐1324‐9.

Pond, S. L. K., A. F. Y. Poon, R. Velazquez, et al. 2020. “HyPhy 2.5‐A Customizable Platform for Evolutionary Hypothesis Testing Using Phylogenies.” Molecular Biology and Evolution 37, no. 1: 295–299. https://doi.org/10.1093/molbev/msz197.

Quinlan, A. R., and I. M. Hall. 2010. “BEDTools: A Flexible Suite of Utilities for Comparing Genomic Features.” Bioinformatics 26, no. 6: 841–842. https://doi.org/10.1093/bioinformatics/btq033.

Ranganathan, Y., M. Ghara, and R. M. Borges. 2010. “Temporal Associations in Fig‐Wasp‐Ant Interactions: Diel and Phenological Patterns.” Entomologia Experimentalis et Applicata 137, no. 1: 50–61. https://doi.org/10.1111/j.1570‐7458.2010.01038.x.

Revell, L. J. 2024. “Phytools 2.0: An Updated R Ecosystem for Phylogenetic Comparative Methods (And Other Things).” PeerJ 12: e16505. https://doi.org/10.7717/peerj.16505.

Roquis, D., M. Robertson, L. Yu, M. Thieme, M. Julkowska, and E. Bucher. 2021. “Genomic Impact of Stress‐Induced Transposable Element Mobility in Arabidopsis.” Nucleic Acids Research 49, no. 18: 10431–10447. https://doi.org/10.1093/nar/gkab828.

Schrader, L., J. W. Kim, D. Ence, et al. 2014. “Transposable Element Islands Facilitate Adaptation to Novel Environments in an Invasive Species.” Nature Communications 5: 5495. https://doi.org/10.1038/ncomms6495.

Sessegolo, C., N. Burlet, and A. Haudry. 2016. “Strong Phylogenetic Inertia on Genome Size and Transposable Element Content Among 26 Species of Flies.” Biology Letters 12, no. 8: 20160407. https://doi.org/10.1098/rsbl.2016.0407.

Stritt, C., S. P. Gordon, T. Wicker, J. P. Vogel, and A. C. Roulin. 2018. “Recent Activity in Expanding Populations and Purifying Selection Have Shaped Transposable Element Landscapes Across Natural Accessions of the Mediterranean Grass Brachypodium distachyon.” Genome Biology and Evolution 10, no. 1: 304–318. https://doi.org/10.1093/gbe/evx276.

Szitenberg, A., S. Cha, C. H. Opperman, D. M. Bird, M. L. Blaxter, and D. H. Lunt. 2016. “Genetic Drift, Not Life History or RNAi, Determine Long‐Term Evolution of Transposable Elements.” Genome Biology and Evolution 8, no. 9: 2964–2978. https://doi.org/10.1093/gbe/evw208.

Talavera, G., and J. Castresana. 2007. “Improvement of Phylogenies After Removing Divergent and Ambiguously Aligned Blocks From Protein Sequence Alignments.” Systematic Biology 56, no. 4: 564–577. https://doi.org/10.1080/10635150701472164.

Venner, S., V. Miele, C. Terzian, et al. 2017. “Ecological Networks to Unravel the Routes to Horizontal Transposon Transfers.” PLoS Biology 15, no. 2: e2001536. https://doi.org/10.1371/journal.pbio.2001536.

Walker, B. J., T. Abeel, T. Shea, et al. 2014. “Pilon: An Integrated Tool for Comprehensive Microbial Variant Detection and Genome Assembly Improvement.” PLoS One 9, no. 11: e112963. https://doi.org/10.1371/journal.pone.0112963.

Walsh, A. M., R. D. Kortschak, M. G. Gardner, T. Bertozzi, and D. L. Adelson. 2013. “Widespread Horizontal Transfer of Retrotransposons.” Proceedings of the National Academy of Sciences of the United States of America 110, no. 3: 1012–1016. https://doi.org/10.1073/pnas.1205856110.

Wang, R. W., B. F. Sun, and Q. Zheng. 2010. “Diffusive Coevolution and Mutualism Maintenance Mechanisms in a Fig‐Fig Wasp System.” Ecology 91, no. 5: 1308–1316. https://doi.org/10.1890/09‐1446.1.

Wang, Y., L. Jin, G. Belusic, L. W. Beukeboom, B. Wertheim, and R. A. Hut. 2024. “Circadian Entrainment to Red‐Light Zeitgebers and Action Spectrum for Entrainment in the Jewel Wasp Nasonia Vitripennis.” Journal of Comparative Physiology A 210, no. 3: 459–472. https://doi.org/10.1007/s00359‐023‐01672‐4.

Wu, C., and J. Lu. 2019. “Diversification of Transposable Elements in Arthropods and Its Impact on Genome Evolution.” Genes (Basel) 10, no. 5: 338. https://doi.org/10.3390/genes10050338.

Xiao, J. H., X. Q. Wei, Y. Zhou, et al. 2021. “Genomes of 12 Fig Wasps Provide Insights Into the Adaptation of Pollinators to Fig Syconia.” Journal of Genetics and Genomics 48, no. 3: 225–236. https://doi.org/10.1016/j.jgg.2021.02.010.

Xiao, J. H., Z. Yue, L. Y. Jia, et al. 2013. “Obligate Mutualism Within a Host Drives the Extreme Specialization of a Fig Wasp Genome.” Genome Biology 14, no. 12: R141. https://doi.org/10.1186/gb‐2013‐14‐12‐r141.

Xie, H. X., X. X. Liang, Z. Q. Chen, et al. 2022. “Ancient Demographics Determine the Effectiveness of Genetic Purging in Endangered Lizards.” Molecular Biology and Evolution 39, no. 1: msab359. https://doi.org/10.1093/molbev/msab359.

Yang, J. S., Z. H. Qian, T. Shi, Z. Z. Li, and J. M. Chen. 2022. “Chromosome‐Level Genome Assembly of the Aquatic Plant Nymphoides Indica Reveals Transposable Element Bursts and NBS‐LRR Gene Family Expansion Shedding Light on Its Invasiveness.” DNA Research 29, no. 4: 1–11. https://doi.org/10.1093/dnares/dsac022.

Yang, Z. H. 2007. “PAML 4: Phylogenetic Analysis by Maximum Likelihood.” Molecular Biology and Evolution 24, no. 8: 1586–1591. https://doi.org/10.1093/molbev/msm088.