Differential LTR-retrotransposon dynamics across polyploidization, speciation, domestication and improvement of cotton (Gossypium)

Differential LTR-retrotransposon dynamics across polyploidization, speciation, domestication and improvement of cotton (Gossypium) | Research Square window.SnipcartSettings = { analytics: { enabled: false } }; (function() { var accessVector = localStorage.getItem('access_vector') || ''; window.dataLayer = window.dataLayer || []; if (accessVector) { window.dataLayer.push({ user: { profile: { profileInfo: { snid: accessVector } } } }); } })(); (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],j=d.createElement(s),dl=l!='dataLayer'?'&l='+l:'';j.async=true;j.src='https://www.googletagmanager.com/gtm.js?id='+i+dl;f.parentNode.insertBefore(j,f);})(window,document,'script','dataLayer','GTM-K279D39R'); Browse Preprints In Review Journals COVID-19 Preprints AJE Video Bytes Research Tools Research Promotion AJE Professional Editing AJE Rubriq About Preprint Platform In Review Editorial Policies Our Team Advisory Board Help Center Sign In Submit a Preprint Cite Share Download PDF Research Article Differential LTR-retrotransposon dynamics across polyploidization, speciation, domestication and improvement of cotton (Gossypium) Lucía Campos-Dominguez, Raúl Castanera, Corrinne E. Grover, Jonathan F. Wendel, and 1 more This is a preprint; it has not been peer reviewed by a journal. https://doi.org/ 10.21203/rs.3.rs-6172192/v1 This work is licensed under a CC BY 4.0 License Status: Published Journal Publication published 27 Oct, 2025 Read the published version in Genome Biology → Version 1 posted 15 You are reading this latest preprint version Abstract Background Transposable elements (TEs) are major components of plant genomes and major drivers of plant genome evolution. The cotton genus ( Gossypium ) is an excellent evolutionary model for polyploidization, speciation, domestication and crop improvement. Here, we implement genome and pangenome analyses to study in detail the dynamics of LTR-retrotransposons (LTR-RT) during the cotton evolution. Results We show that some LTR-RT lineages amplified in tetraploid cotton compared to their diploid progenitors, whereas others stayed stable or amplified but were removed through solo-LTR formation. Using species-level pangenomes we show that only a few lineages (CRM, Tekay, Ivana and Tork) remained active after polyploidization and are still transposing. Tekay and CRM elements have re-shaped the centromeric and pericentromeric regions of tetraploid cottons in a subgenome specific manner, through new insertions but also selective eliminations through solo-LTR formation. On the other hand, Ivana and Tork have actively inserted within or close to genes. Finally, population-level analyses using the two pangenomes and data from 283 and 223 varieties of G. hirsutum and G. barbandense reveal changes in Transposon Insertion Polymorphism (TIP) frequencies accompanying domestication and improvement of both species, suggesting the possibility of selection on linked regions. Conclusions Our findings reveal that LTR-RT lineages followed differential dynamics during cotton evolution, displaying differences among species and the two coresident genomes of allopolyploid cotton. A handful of the LTR-RT lineages that expanded after polyploidisation helped shape the genomes of both G. hirstutum and G. barbadense , impacting their centromere and pericentromeric regions as well as protein- coding genes. retrotransposon solo-LTR allopolyploid plants pangenome TIP (Transposon Insertion Polymorphism) evolutionary genomics Figures Figure 1 Figure 2 Figure 3 Figure 4 Figure 5 Figure 6 Background Transposable elements (TEs) are major components of plant genomes. Their repetitive nature and their ability to amplify and move within the genome makes them a major source of genetic variability and one of the major drivers of plant genome evolution [ 1 , 2 ]. TEs are an important source of mutations underlying crop domestication and improvement [ 3 ]. However, as potential mutagens TEs are tightly controlled, mainly through epigenetic mechanisms [ 4 ]. For this reason, the dynamics of TEs within genomes has been often seen as an arms-race between an invasive parasite and a defensive host, with genomes tightly controlling TE activity and TEs evolving to evade genomic suppression [ 5 – 7 ]. TE activity is known to increase in situations where DNA methylation and epigenetic silencing decrease, such as stress situations [ 8 ], including the genomic shock that results from the combination of two different genomes in polyploid species [ 9 , 10 ]. This activation can lead to a burst of transposition and a subsequent repression when silencing is re-established [ 11 ]. Yet accumulated data during the last few years in plant genomes suggests that the dynamics of different TE types, lineages and families that coexist in a genome are highly diverse. These lineages and families may thus be viewed as different players in an evolving ecosystem [ 12 – 14 ]. Long terminal repeat retrotransposons (LTR-RTs) are the most prevalent type of TEs in plants, as exemplified by maize, where LTR-RTs account for 90% of the TEs and more than 80% of the genome space [ 14 ]. LTR-RTs are divided into two major superfamilies, Copia and Gypsy , which often show important differences in their regulation and distribution within genomes [ 15 ]. LTR-RTs transpose through a replicative mechanism, i.e., via an RNA copy reverse transcribed into cDNA prior to integration at a new location. Therefore, their activity results in an increase in copy number, which endows them with the potential to invade, and thereby inflate, the genome. As an example, the genome of the wild rice O. australiensis doubled in size in just three million years by the amplification of only three families of LTR-RTs [ 16 ]. Counterbalancing this genome size expansion via TE proliferation, LTR-RTs can also be eliminated from genomes by different mechanisms [ 17 , 18 ]. The main mechanism for this deletional process is illegitimate recombination, often at the LTRs, giving rise to so-called solo-LTRs resulting in truncated LTR-copies [ 19 – 22 ]. Although the importance of recombination-mediated LTR-RT elimination is well recognized, its analysis is often neglected when studying LTR-RT dynamics. The cotton genus ( Gossypium ) serves as an excellent evolutionary model for polyploidization, domestication and crop improvement. Allopolyploid cotton species arose 1–2 mya from the interspecific hybridisation of an “A-genome” and a “D-genome” species, which diverged from each other 5–7 mya [ 23 ]. Because of this long divergence of the two genome donors, A and D, their resulting subgenomes can be easily distinguished in the tetraploid genomes. The ancestor providing the A subgenome is closely related to the modern species G . herbaceum (A1) and G. arboreum (A2), whereas the ancestor of the D subgenome is closest, among extant D-genome species, to G. raimondii (D5) [ 23 ]. Two allotetraploid species ( G. hirsutum , AD1, and G. barbadense , AD2) were independently domesticated from different wild progenitors in Central and South America, respectively, over the last 8000 years [ 23 , 24 ]. TE studies on cotton species have previously highlighted their significant role on genome size expansion, being the major contributors to genome size differences between diploid cotton species [ 25 ]. Recent analyses have also suggested that TEs may have amplified accompanying the polyploidization event, particularly for the Gypsy Tekay and CRM lineages of LTR-RTs [ 26 ], which may have reshaped the centromeres in tetraploids [ 27 ]. Here we study the LTR-RT dynamics (both LTR-RT insertion and removal through solo-LTR formation) accompanying diploid divergence, polyploidization, and speciation at the allopolyploid level using genome assemblies of eight G. hirsutum and ten G. barbadense accessions and representatives of their parental diploid species. Using pangenomes generated for the two domesticated tetraploid cotton species ( G. hirsutum and G. barbadense ), we report on recent LTR-RT dynamics accompanying initial crop domestication and subsequent improvement. Using this temporally stratified and explicit approach, we show that each LTR-RT lineage has experienced its own dynamics throughout the different events characterizing the recent evolution of cotton species and varieties. We also use the pangenomes to analyze relatively recent or ongoing transposition events within the two allopolyploid species, documenting the frequency of the intraspecies transposon insertion polymorphisms (TIPs) in approximately 500 sequenced genomes representing the wild-to-domesticated continuum within G. hirsutum and G. barbadense. Our results show that LTR-RTs may have been instrumental in shaping cotton genomes and generating the genome variability that was subject to selection during domestication and crop improvement. Results LTR-RT lineages have differentially amplified in tetraploid cotton species As a first step in exploring the dynamics of LTR-RTs in the recent evolution of cotton, we annotated intact LTR-RT copies in the highly contiguous assemblies of eight accessions of G. hirsutum and ten accessions of G. barbadense , as well as in the genomes of two diploid species that represent the best models of the ancestral A and D subgenomes of cotton ( G. herbaceum and G. raimondii , respectively). Our analysis focused only on the intact LTR-RT copies (identified by stringent approaches; see methods) to concentrate on the youngest elements, i.e., those that could have arisen from transposition accompanying polyploidization, speciation, and domestication of cotton, events that occurred within the last 1–2 million years. Results in Fig. 1 a (left) show that the two polyploid species have, in general, a higher number of intact LTR-RTs per megabase than do the diploids, especially in the D subgenome, where the density of intact LTR-RT elements is roughly doubled (from 3.5 elements per Mbp in the D-genome diploid parental G. raimondii to ~ 8 and ~ 7 elements per Mbp in G. hirsutum and G. barbadense , respectively). For the A subgenome, there is an increased density with respect to that of G. herbaceum but not with respect to G. arboreum . However, it has been shown that G. arboreum experienced a recent increase in LTR-RT copy number, and therefore its present LTR-RT density does not reflect its density when polyploidization occurred [ 28 ]. Therefore, we only use G. herbaceum as a comparator for the A subgenomes henceforth. Between the two tetraploid species, G. hirsutum varieties generally have a slightly higher number of intact LTR-RTs than the G. barbadense accessions, and there are no major differences among accessions within each of the two species (Fig. 1 a, left panel). Insertion times were calculated for each element by comparing the two LTR sequences for each intact LTR-RT. Figure 1 a (right) shows that LTR-RT insertions are younger, on average, in the polyploids than in the diploids (Fig. 1 a, right panel), suggesting more recent activation of LTR-RTs in the polyploids. Analysis of the LTR-RT composition of the genome of the two tetraploid species relative to their two diploid parentals reveals that the increase in intact copies, and therefore in the percentage of the genome covered by LTR-RTs, is mainly due to few LTR-RT lineages (e.g., CRM and Tork, see Fig. 1 b). Analysis of the number of intact elements from each LTR-RT lineage (Figs. 2 , 3 and Supplementary Fig. 1) indicates that there is diversity in dynamics even among lineages belonging to the same superfamily of LTR-RTs (i.e. Copia and Gypsy ), whereby some lineages exhibit abundant increases in copy number in the polyploids whereas others remain relatively stable during recent cotton evolution. For example, copy number for the Copia Tork lineage is three times higher per megabase in the two tetraploids than in the two diploid progenitors (Fig. 2 top left). Insertion time analysis of Tork elements shows that although most elements are relatively old in the diploids, they are notably younger in the polyploids (Fig. 2 middle left). This confirms that Tork has been particularly active in the tetraploid genomes and may still be actively transposing. In contrast, the number of copies per megabase of other lineages, such as the relatively abundant Gypsy Athila (Fig. 2 top right), is similar when comparing each subgenome of the polyploid with its corresponding diploid genome, suggesting no or little active retrotransposition since polyploid formation. An analysis of the insertion time of Athila intact elements shows a relatively old insertion time, particularly for elements residing in the A subgenomes, and an almost identical distribution of insertion times between diploid and tetraploid genomes, supporting very low activity after polyploidization (Fig. 2 middle right). However, Athila elements have probably been very active in the past, as Athila-related sequences account for a sizeable fraction (~ 13% and 5% in subgenomes A and D, respectively) of the genome of both tetraploid species and their diploid parentals (~ 10% in G. herbaceum and 4,6% in G. raimondii ) (Fig. 1 b). A general inference from the analysis of copy number for all LTR-RT lineages is that the evolutionary behavior of the different TE families varies considerably among lineages, subgenomes, species and, with respect to timing of proliferation events (Fig. 2 , 3 and Supplementary Fig. 1). Some LTR-RT lineages show an increase in retrotransposition in both subgenomes after polyploidisation (the Copias Tork and Ivana, and the Gypsy CRM), whereas others mainly increased in one subgenome (the Gypsy Tekay). In contrast, a stable number of LTR-RT insertions was observed among all species at both the diploid and tetraploid level for some retrotransposon lineages, such as the Gypsy Athila, Ogre, Reina, Galadriel families and the Copia Ale family. Finally, some TE families are present at similarly low copy numbers in all genomes (the Copias Sire, Bianca, Ikeros and TAR). In general, lineages present at a relatively high copy number exhibited similar insertion numbers when comparing different varieties of the same species, although this is not always the case, as evidenced by the Ivana lineage in the subgenome D of G. hirsutum (Supplementary Fig. 1). Similarly, the copy number of each lineage is generally similar between the two tetraploid species, although differences for specific lineages do exist. Notably among these exceptions are the CRM and Tekay Gypsy lineages (Fig. 3 ), which are further discussed below. The rate of LTR-RT elimination varies among lineages, species and subgenomes We further explored the dynamics of LTR-RT elements accompanying recent cotton evolution by quantifying LTR-RT elimination though intra element LTR-recombination. We annotated potential solo-LTRs by searching for LTRs corresponding to the lineages and families of intact LTR-RTs that did not have surrounding sequences that could correspond to LTR-RT internal sequences (see methods). These results show that there are important differences among the LTR-RT lineages (Fig. 2 , 3 and supplementary Fig. 1 bottom panels). For example, the relative number of Tork solo-LTRs is similar between tetraploid and diploid genomes (Fig. 2 bottom left), despite the increase in intact Tork insertions in polyploid cotton relative to their parental diploids (Fig. 2 top left). Combined with the younger age of intact Tork elements in the tetraploids (relative to the diploids; Fig. 2 middle left), these observations suggest both that Tork elements have been more active in the tetraploids than in the diploids (after polyploidization) and that only a minor fraction of the new insertions has been removed by recombination. In contrast, both diploids and tetraploids contain two to ten times more Athila solo-LTRs than intact elements (Fig. 2 right), suggesting a high rate of elimination of these elements. As noted above, Athila has similarly old element insertion times between diploids and polyploids genomes (Fig. 2 middle right), suggesting low or no recent activity of this family; however, the greater number of solo-LTRs in G. barbadense relative to G. hirsutum may indicate recent activity and rapid turnover of Athila insertions in G. barbadense . A similar diversity of solo-LTR dynamics characterizes the remaining LTR-RT lineages. For example, Ogre solo-LTRs are also more abundant than intact elements (supplementary Fig. 1); however, unlike Athila, intact elements tend to be younger in the tetraploids, particularly for those in the D subgenome of G. hirsutum , where the number of solo-LTRs is particularly high. This observation may suggest that the removal of intact elements through solo-LTR formation affected mostly older elements. On the other hand, Ivana, very much like Tork, shows a low number of solo-LTRs and a high number of young intact elements, suggesting active retrotransposition and low recombination rates. In contrast, Angela, Sire, TAR, Bianca and Ikeros families have more solo-LTRs than intact elements, although the numbers are all low, suggesting limited retrotransposition activity for each of these families. Finally, Ale, Galadriel, and Reina show the same number of intact LTR-RTs and solo-LTRs in the tetraploids as in the parental diploids, and similar LTR-RT age, suggesting little or no activity during the recent evolution of Gossypium . CRM and Tekay insertion and elimination have shaped the centromeres of tetraploid cottons The CRM and Tekay lineages have previously been shown to have amplified in the tetraploids [ 26 ] and have been suggested to play a role in reorganization of their centromeres [ 27 ]. Our results confirm that CRM and Tekay have been actively transposing after the polyploidization event, as deduced by their recent insertion time distribution (Fig. 3 ). Moreover, our results also show that CRM copy number increased greatly in the A subgenomes of both tetraploids (up to 8-fold in G. hirsutum and up to 6.5-fold in G. barbadense ) and to a lesser extent in the D subgenomes (up to 2.6-fold in G. hirsutum and up to 2-fold in G. barbadense ), with this increase being slightly higher in G. hirsutum than in G. barbadense (Fig. 3 top left). A different pattern was observed for Tekay elements, which exhibited increased copy numbers in the D subgenomes of both tetraploids (also slightly more in G. hirsutum than in G. barbadense ), but not in the A subgenomes, where its copy number remained high and stable in G. hirsutum and decreased about 30% in G. barbadense (Fig. 3 top right). These results thus show that, whereas both CRM and Tekay actively transposed accompanying polyploidization, their amplification rate was dissimilar and differed between species and subgenomes. Analysis of the CRM and Tekay solo-LTRs shows that the elimination rate also varies between the two polyploids and between subgenomes (Fig. 3 bottom panel). The number of CRM solo-LTRs is much higher in the D subgenomes than in the A subgenomes in both polyploids, but it is particularly high in G. barbadense . This could be the consequence of a higher elimination rate of recent CRM insertions from the D subgenomes, particularly in G. barbadense , which would also explain the lower number of intact CRMs in the D subgenomes of this species (Fig. 3 top left) and the slightly older insertion time of the remaining G. barbadense CRMs (Fig. 3 middle panel). Interestingly, the analysis of the Tekay lineage suggests an opposite dynamic. Tekay solo-LTRs are much more abundant in the A subgenomes of both tetraploids than in the D subgenomes, a trend that is even more extreme for G. hirsutum (Fig. 3 bottom right). Interestingly, the high number of Tekay solo-LTRs, doubling that of complete elements, is also seen in the diploid parental, G. herbaceum , which may suggest a link between the nature of the A subgenomes and the efficiency of Tekay recombination. This high elimination rate of Tekay elements in the A subgenomes could explain why the recent activity suggested by the very recent age of the intact elements (Fig. 3 middle panel) has only resulted in an increase of Tekay insertion in the D subgenomes (Fig. 3 top panel). The distribution of the intact CRM elements and solo-LTR CRM shows a high concentration of CRM intact copies in the regions characterized as containing the centromeres [ 27 ], whereas CRM solo-LTRs are more homogeneously distributed along the chromosomes (Fig. 4 a, Supplementary Fig. 2), suggesting that, although CRM elements may target the centromere regions [ 29 ], their differential elimination from other chromosome regions reinforces their concentration in the centromeric regions. Similarly, intact Tekay elements in the D subgenomes are concentrated in the pericentromeric regions whereas Tekay solo-LTRs are not, suggesting active elimination of Tekay elements from chromosomal regions maintains their prevalence in pericentromeric regions. Interestingly, Tekay elements are not found in CRM-rich regions of the centromere, perhaps reflecting an insertional preference difference between these two types of elements. The distribution of Tekay intact and solo-LTRs in the A subgenomes seems less skewed than in the D subgenome, although they also show an opposite distribution with respect to genes. Comparison of the CRM and Tekay distributions between the diploids and the tetraploids also reveals notable differences. First, the high increase in CRM copy number in the A subgenomes of the two polyploids, and particularly in G. hirsutum (Supplementary Fig. 3), is concentrated almost exclusively in the centromere, whereas in the diploids, the distribution of the CRM elements is less skewed (Fig. 4 b, and Supplementary Fig. 3). Second, the solo-LTR distribution in the diploids parallels that of the intact copies, suggesting that there is no preferential elimination of CRM elements from chromosome arms (Fig. 4 b). For Tekay elements, the increase in intact elements in the D subgenomes of polyploids (up to 5-fold in G. hirsutum and up to 3.5-fold in G. barbadense) without a concomitant increase in number of solo-LTRs appears to concentrate the intact elements in the pericentromeric regions (Fig. 4 b). A small number of LTR-RT lineages accounts for intraspecies polymorphisms Although the number of insertions of most LTR-RT lineages appears relatively constant within each tetraploid species, the young age of the insertions for many LTR-RT lineages suggests recent retrotransposition. Therefore, to understand the scale and scope of TE activity within species, we constructed reference-based pangenomes for both species that we subsequently used to characterize transposon insertion polymorphisms (TIPs) within G. hirsutum and G. barbadense . We used genome assemblies of seven different varieties to build the G.hirsutum pangenome and genome assemblies of nine varieties and nine different varieties to build the G. barbadense pangenome. The different genomes sequences were aligned to the references for each species (TM-1 and 3–79) to build each pangenome. The alignments covered 94.1–98.7% the G. hirsutum TM1 reference and 94.8–97.0% of the G. barbadense 3–79 reference. The number of structural variants (SV) detected in the G. hirsutum and G. barbadense pangenomes was 63,081 and 55,122, respectively, most of which (57,727 and 48,990, respectively) were insertion/deletion (indel) polymorphisms. Approximately 20% of these indels (12,316 in G. hirsutum and 9,386 in G. barbadense ) contained LTR-RT sequences (at least 30% coverage of an LTR-RT). To detect recent LTR-RT insertions, we defined LTR-RTs TIPs as indel SVs with at least 80% identity over 80% of the length of an intact LTR-RT representative or a solo-LTR, and a minimum size of 80% of an intact LTR-RT or a solo-LTR. This analysis led to the detection of 4436 LTR-RT TIPs in G. hirsutum and 3520 in G. barbadense , and 2406 solo-LTR TIPs in G. hirsutum and 1313 in G. barbadense . These SVs were homogeneously distributed along chromosomes, as were the LTR-RTs and solo-LTR TIPs (Supplementary Fig. 4a). In addition to this overall quantitative view of TE-mediated TIPS, LTR-RT lineage classification of both intact and solo-LTR TIPs shows that only few of the many LTR-RT lineages are responsible for the polymorphism generated within the two cotton species (Figs. 5 a and 5 b, Supplementary Tables 1,2). Interestingly, some lineages that experienced copy number increase in tetraploid cotton relative to their diploid progenitors, e.g., Ogre (Supplementary Fig. 1), also show little variability within the polyploids (less than 1%, Supplementary Tables 1,2 and Fig. 5 a), suggesting that they experienced a burst of transposition following polyploidization but became silent before the split of the two tetraploids from their common ancestor. In contrast, some lineages appear to have retained high transpositional activity after speciation of the two polyploids. These include the Gypsy CRM and Tekay lineages and the Copia Ivana and Tork lineages (Fig. 5 a). Interestingly, CRM is mainly polymorphic as solo-LTRs, whereas Ivana and Tork are mainly polymorphic as complete LTR-RTs and Tekay is notably polymorphic as both. These patterns could suggest recent elimination activity for CRM, recent retrotransposition for Ivana and Tork, and insertion followed by elimination of new copies for Tekay. Whereas CRM and Tekay elements concentrate in centromeric and pericentromeric regions (see above), Copia elements Ivana and Tork have a more widespread chromosomal distribution (not shown) and are more likely to have been inserted near genes and regulatory regions. We explored the potential effect of these polymorphic Tork and Ivana insertions by evaluating their position relative to genes in the pangenomes of both species (Fig. 5 b). Our analysis shows that Ivana and Tork insertions are significantly enriched inside and near (2 kbp) genes in both species and may therefore have impacted gene coding capacity and regulation, potentially triggering phenotypic consequences. LTR-RT TIPs impact on cotton domestication and breeding To address the possibility that genomic variability triggered by recent TIPs had an impact on cotton domestication and improvement, we analyzed their presence in wild and domesticated G. hirsutum and G. barbadense populations, including landraces and elite varieties. To this end we mapped short-read resequencing data of 283 varieties of G. hirsutum and 223 varieties of G. barbadense to their respective pangenome graphs and genotyped the presence/absence of the 6841 G. hirsutum and 5749 G. barbadense TIPs across these populations. The varieties analyzed represent the wild-to-domesticated continuum of both species including, 151 and 51 elite varieties, 66 and 60 landraces, and 73 and 112 wild accessions, of G. hirsutum and G. barbadense , respectively (Supplementary Table 3) [ 24 ]. Principal Component Analysis (PCA) shows that TIPs provide enough signal to distinguish between groups of wild, landrace, and domesticated varieties and suggests that there has been extensive differential transposition activity and/or retention/elimination of LTR-RT copies and solo-LTRs between the three groups of accessions (Fig. 6 a,b). This is particularly clear for G. hirsutum , in line with previous analyses based on SNPs which showed that the cultivated G. hirsutum accessions clustered tightly due to their narrow genetic diversity [ 24 ]. This drastic reduction of diversity accompanying domestication was not seen in G. barbadense , probably because modern G. barbadense cultivars have a more complex and obscure origin that in addition involves many intentional introgressions [ 24 ]. An analysis of the TIP frequency in G. hirsutum and G. barbadense populations shows that there are few LTR-RT insertions or deletions present in only one or two groups (see Supplementary Table 4) whereas most TIPs (5860 in G. hirsutum and 4450 in G. barbadense ) are present in all three groups (i.e., wild, landrace, and domesticated). Interestingly, however, the TIPs shared by the three groups are present in very different frequencies among accession in those groups (Supplementary Figs. 5 and 6). The presence of TIPs at different population frequencies suggests the possibility of selection, positive or purifying, acting on polymorphic LTR-RTs or on closely linked loci . To find signs indicative of positive selection of TIPs in cultivated varieties, we used the population branch statistics (PBS) method. This approach measures the level of genetic differentiation for a specific population along a phylogenetic branch, providing insights into population-specific evolutionary changes [ 30 ]. Figure 6 c and 6 d show the PBS distribution of LTR-RT insertions and deletions relative to the reference genome in cultivated versus landrace and wild accessions and in wild versus cultivated and landrace accessions in G. hirsutum. We identified the TIPs with the highest PBS (above 95th percentile in the comparisons cultivated against landraces and wild, and wild against cultivated and landraces). Eighty percent of the high-PBS TIPs of the wild versus cultivated and landraces (linked to domestication) and 63% of high-PBS TIPs in cultivated versus wild and landraces (linked to cotton breeding) are related to Gypsy LTR-RT (mainly Tekay) elements (Supplementary Table 5). Tekay elements are concentrated in pericentromeric region of chromosomes, particularly in the D subgenomes, and our results show that they have contributed to the recent evolution of these regions, as discussed above. Many of the Tekay TIPs with high PBS are far from genes (Supplementary Table 4), have a distribution similar to the Tekay intact elements (Fig. 4 and Supplementary Fig. 4b), and are in regions characterized by a low density in coding genes. The results presented here suggest that Tekay insertions in these regions may have been selected in the recent evolution of G. hirsutum and G. barbadense . Notably, although most high-PBS Tekay TIPs are far from genes, others are located close to genes (less than 2 Kb away, Supplementary Table 5 and Suplementary Fig. 4 b) and may have impacted their coding capacity or regulation. This is also the case of High-PBS TIPs related to Copia LTR-RTs, mainly the Ivana and Tork superfamilies, which are frequently found tightly associated with genes. The high PBS associated with these TIPs may suggest selection due to favorable or detrimental effects on the nearby genes. We analyzed the expression of genes tightly linked with high-PBS TIPs using previously published data for the G. hirsutum TM1 accession [ 31 ] and found that they show different patterns of expression throughout cotton development and under different stress conditions (Supplementary Fig. 7). An analysis of the functions of these genes shows a significant enrichment for "auxin response" genes (p-value = 0,00039, GO term GO:0009733) in G. hirsutum . We found 5 genes belonging to this functional category (four coding for a SAUR-like auxin-responsive proteins and one for an Indoleacetic acid-induced protein 16) within the 95 percentile of highest PBS in G. hirsutum , indicative of high population frequency differentiation. This represents an 11-fold enrichment in comparison to the total annotated predicted genes (Fig. 6 e). Two of these TIPs are shown in Fig. 6 g and 6 h. For the first of these, we observed a strong population frequency differentiation of the landraces and cultivated accessions with respect to the wild (potential relationship with domestication), whereas for the second, we found a strong differentiation between the cultivated and the other two groups, suggesting that this TIP could have been targeted during crop improvement. A preliminary analysis of the TIPs with high PBS value in G. barbadense TIPs also showed an enrichment for the function "response to auxin” (p = 0.022, Fig. 6 e, Supplementary Table 4), suggesting that variants in auxin-related genes were selected in the domestication and breeding process of the two species of cotton. Whether this result reflects parallel selection or historical interspecific introgression is at present an open question. Discussion Plant TEs, and in particular LTR-RTs, have been traditionally assumed to evolve through bursts of amplification followed by periods of low activity and elimination [ 11 ]. Recent data, however, suggest that the diversity of LTR-RTs that coexist within plant genomes often results from heterogeneous evolutionary dynamics [ 12 ]. The results presented here show that LTR-RT elements are in general more abundant and recent in the two tetraploid cotton species compared with models of their diploid progenitors, suggesting a general burst of amplification accompanying the polyploidization event, as it is often but not always the case in plants [ 9 ]. However, these results also show that different LTR-RT lineages experienced rather different evolutionary dynamics during cotton polyploidization, speciation, and subsequent crop improvement processes. As an example, whereas Ogre elements amplified in the tetraploids, they probably became silent before speciation, as they show little polymorphism between the two tetraploid species. Other LTR-RT lineages, such as Tork and Ivana, have probably retained activity for longer, since more than half of their insertions are polymorphic within the species. This is also the case for CRM and Tekay, which have recently been proposed to play important roles in the evolution of Gossypium centromeres [ 27 , 28 ]. Our results show that CRM elements amplified after polyploidization. This CRM element number increase has been higher in the A than the D subgenomes, probably because they are more efficiently removed from the former subgenome, as evidenced by a higher number of solo-LTRs. As CRMs concentrate in cotton centromeric regions, this leads to a differential increase of CRMs in the centromeric regions of the A subgenomes in the two tetraploid species as compared to the model diploid progenitor genome (represented by modern G. herbaceum ). Moreover, the number of CRM solo-LTRs in the D subgenomes is higher in G. barbadense than in G. hirsutum , which correlates with a lower number of intact CRMs in this species, suggesting that CRM elimination from the centromeric regions of the D subgenomes is higher in this species. In contrast to CRM LTR-RTs, Tekay elements have greatly increased in the D subgenomes while their number has remained constant in the A subgenomes. These results are in line with recent evidence regarding the role of these two lineages in centromere evolution in cotton after polyploidization, with the proportion of CRM versus Tekay elements being higher in the centromeric regions of the A subgenomes, and lower in those of the D subgenomes [ 27 ]. Our results add to these earlier observations in suggesting that this shift may reflect differences in elimination of the insertions between the two subgenomes rather than differential insertion. Moreover, the analysis of CRM and Tekay intact and solo-LTR distributions along chromosomes suggests that, although these elements may be centrophilic[ 32 ] and target the centromere for integration, their elimination via solo-LTR formation from chromosome arms reinforces their concentration in centromeric and pericentromeric regions. Taken together, our results are compatible with the important role of CRM and Tekay LTR-RTs in the expansion and reorganization of cotton centromeres after polyploidization that has been recently suggested[ 27 , 28 ] and points to the role of their selective elimination from chromosomal regions and subgenomes in specifying their dynamics and their genomic distribution. With respect to continuing LTR-RT activity within species, our results show that CRM and especially Tekay intact LTR-RT and solo-LTR insertions are highly polymorphic within the two polyploid species, and are present at different population frequencies in wild, landrace and cultivar groups of accessions. Interestingly, many of the Tekay TIPs are among those with the highest Population Branch Statistic (PBS) values, suggesting potential selection during domestication and breeding. The potential selection of Tekay TIPs sitting in pericentromeric regions suggests that variants of these heterochromatic regions could have selectable phenotypic impacts. Centromeres are fast evolving structures[ 29 ] and it has previously been shown that centromere variants have been selected during the domestication of crops, such as maize and wheat [ 33 , 34 ]. The selection of centromere variants could be both due to differences in the centromere function or due to the selection of particular alleles of genes residing in this low recombining region. Clarifying this point will require further analyses, but the results presented here point to an important role of genetic variability in heterochromatic regions for crop domestication and breeding. In addition to being important players in centromere evolution [ 29 ], LTR-RTs have also played a major role as drivers of genetic diversity in plant gene coding capacity and regulation [ 2 , 35 ]. Our data show that, in addition to Gypsy CRM and Tekay, two Copia LTR-RT lineages (Tork and Ivana) have amplified and remained active in the two polyploid species of cotton since their divergence from a common ancestor. Notably, we identified a large number of very recent retrotranspositions, as evidenced by their identical LTRs, that are significantly enriched within exons and in the proximal regions upstream and downstream of genes. Many of these insertions have varying frequencies in the wild, landrace, and cultivar groups, exhibiting high PBS values and suggesting that selection during cotton domestication and crop improvement may have shaped their distribution. Gene expression data reveals a wide range of developmental and stress-related expression patterns. Considering the many different traits targeted by domestication and crop improvement, including plant architecture, fruiting habit, flower and seed development and fiber-length and quality [ 36 ], more work is needed to determine whether the genomic changes introduced these genic and gene-adjacent LTR-RT insertions have contributed to these processes. It is interesting to note that both G. hirsutum and G. barbadense exhibited an overrepresentation of auxin-responsive genes, and in particular SAUR-like genes, among those impacted by LTR-RT insertions. Auxins are well-established as key controllers of plant development and stress responses [ 37 ], which also play a role in the regulation cotton fiber development [ 38 , 39 ] and broader processes such as plant growth [ 40 , 41 ]. Moreover, SAUR genes, which are key players in plant adaptation and growth [ 42 ] have been shown to have been targets of domestication and breeding of crops such as rice, citrus, and Brassica oleracea [ 43 – 45 ]. The potential independent selection of transposon insertion polymorphisms closely linked to different auxin-related genes in both G. hirsutum and G. barbadense , highlights the importance of this hormone signaling in cotton domestication and breeding and suggests convergent selection in both species (notwithstanding the possibility of an interspecific introgressive origin). Conclusion This study highlights the importance of considering lineage level TE classifications when studying LTR-RT dynamics, as well as examining not only retrotransposition but also recombination and the formation of solo-LTRs as the joint proximal determinants of TE presence in genomes. Our study shows LTR-RT dynamics are highly differential among LTR-RT families, and that only some were activated after polyploidisation and disparately between the two cotton species and the two subgenomes in each. CRM and Tekay elements appear to have played key roles in centromere reorganization, and Tork and Ivana are suggested to actively generate variability within and close to genes. These few recently active lineages have generated LTR-RT polymorphisms that may have been selected during cotton domestication and crop improvement, have modified chromosome heterochromatic regions and have impacted genes expressed in different tissues and environmental situations, including genes related to auxin signaling, which are key players in plant and fiber development regulation, in both G. hirsutum and G. barbadense . Methods All scripts and pipelines used for this study are available within the following directory: https://github.com/Lcamdom/panTEvo Annotation of intact LTR-RTs and solo LTRs All genomes selected and included in this study (supplementary table 6) were assembled into pseudo-chromosomes. These were annotated for intact LTR-RT elements using EDTA-raw scripts on the LTR-mode [ 46 ]. The intact elements found were further classified using TEsorter [ 47 ]. All intact elements within the same species and lineage were clustered using CD-hit [ 48 ] for 80% similarity and 80% coverage. The representatives for each cluster were used to build species-specific, home-made in silico libraries of intact LTR-RT elements. For each element in these libraries, the left LTR was extracted using the structural annotation generated by LTR_Retriever [ 49 ], and species-specific home-made libraries of LTRs were generated. Solo-LTRs were then identified in each genome using these libraries. Blast [ 50 ] was used to find LTRs across the genome, and solo-LTRs were then defined by filtering out any hit with any LTR-RT sequence 1kb either upstream or downstream. Pangenome construction and TIP characterisation Minimap2 [ 51 ] was used to map every G. hirsutum and G. barbadense genome to the references (TM1 and 3–79, respectively). Structural variants (SVs) were called using svim-asm [ 52 ] and the vcfs were merged using bcftools [ 53 ] and collapsed with truvari [ 54 ]. The pangenome graphs were generated using the vg toolkit [ 55 ]. TIPs were characterised using blast [ 50 ] to find homology between the pangenome SVs and the intact LTR-RT and solo-LTR libraries. Resulting hits were excluded when they represented less than 80% of the length of the TE consensus, had lower than 80% homology with the reference TE, and covered less than 80% of the SV. Pangenome genotyping and population frequency analysis Re-sequencing data from 295 and 222 accessions across the wild-to-domesticated continuum (Supplementary Table 3) were mapped to the G. hirsutum and G. barbadense pangenomes, respectively. Mapping of reads to the pangenome graph and structural variant calling were accomplished using the vg toolkit [ 55 ]. The TIP population frequencies were calculated using R. TIPs present in the three population groups were combined in a single matrix to obtain PBS values. TIPs with MAF > 0.01 and call rate > 95% were used to calculate Fst values using SNPready [ 56 ]. PBS values were calculated following the formulas described in [ 30 ]. GO enrichment of genes with proximal TIPs with high PBS was performed with GOATOOLS [ 57 ]. Expression analysis of genes hear high PBS TIPs G. hirsutum RNA-seq data [ 31 ] was processed to estimate transcript abundance using Trinity Transcript Quantification scripts [ 58 ] and Salmon [ 59 ] was used as the abundance estimation method. The raw count matrices were processed into Rlog values using DESeq2 [ 60 ]. Heatmaps are normalised per row (scale = rows), and were generated using the Pheatmap function in R. Declarations Author Contribution LCD and JMC designed the project with the help of all other authors. LCD did most of the analyses with the help of RC for TIP genotyping and TIP population analyses. JFW and CEG contributed cotton sequence information and resequencing data for cotton varieties. LCD and JMC wrote the paper with contributions ofrom all other authors. Data Availability The source code supporting this study is available in the GitHub repository https://github.com/Lcamdom/panTEvo. References Wendel JF, Jackson SA, Meyers BC, Wing RA. Evolution of plant genome architecture. Genome Biol [Internet]. 2016;17:37. Available from: http://genomebiology.biomedcentral.com/articles/10.1186/s13059-016-0908-1 Lisch D. How important are transposons for plant evolution? Nat Rev Genet [Internet]. 2013;14:49–61. Available from: http://www.ncbi.nlm.nih.gov/pubmed/23247435 Andersson L, Purugganan M. Molecular genetic variation of animals and plants under domestication. Proc Natl Acad Sci U S A. 2022;119:e2122150119. Liu P, Cuerda-Gil D, Shahid S, Slotkin R. The Epigenetic Control of the Transposable Element Life Cycle in Plant Genomes and Beyond. Annu Rev Genet. 2022;56:63–87. Luo S, Zhang H, Duan Y, Yao X, Clark AG, Lu J. The evolutionary arms race between transposable elements and piRNAs in Drosophila melanogaster. BMC Evol Biol [Internet]. 2020;20:14. Available from: https://doi.org/10.1186/s12862-020-1580-3 Kosuge M, Ito J, Hamada M. Landscape of evolutionary arms races between transposable elements and KRAB-ZFP family. Sci Rep [Internet]. 2024;14:23358. Available from: https://doi.org/10.1038/s41598-024-73752-7 Lawlor MA, Ellison CE. Evolutionary dynamics between transposable elements and their host genomes: mechanisms of suppression and escape. Curr Opin Genet Dev [Internet]. 2023;82:102092. Available from: https://www.sciencedirect.com/science/article/pii/S0959437X23000722 Zhang H, Lang Z, Zhu JK. Dynamics and function of DNA methylation in plants. Nat Rev Mol Cell Biol. 2018. pp. 489–506. Vicient CM, Casacuberta JM. Impact of transposable elements on polyploid plant genomes. Ann Bot [Internet]. 2017;120:195–207. Available from: http://dx.doi.org/10.1093/aob/mcx078 Nieto Feliner G, Casacuberta J, Wendel JF. Genomics of Evolutionary Novelty in Hybrids and Polyploids. Front Genet. 2020;11:792. Baidouri M, El, Panaud O. Comparative genomic paleontology across plant kingdom reveals the dynamics of TE-driven genome evolution. Genome Biol Evol [Internet]. 2013;5:954–65. Available from: http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3673626&tool=pmcentrez&rendertype=abstract Pulido M, Casacuberta JM. Transposable element evolution in plant genome ecosystems. Curr Opin Plant Biol. 2023;75:102418. Stritt C, Thieme M, Roulin AC. Rare transposable elements challenge the prevailing view of transposition dynamics in plants. Am J Bot. 2021;108:1310–4. Stitzer MC, Anderson SN, Springer NMV, Ross-Ibarra J. The Genomic Ecosystem of Transposable Elements in Maize. PLoS Genet. 2021;17:e1009768. Pereira V. Insertion bias and purifying selection of retrotransposons in the Arabidopsis thaliana genome. Genome Biol. 2004;5:R79. Piegu B, Guyot R, Picault N, Roulin A, Saniyal A, Kim H, et al. Doubling genome size without polyploidization: Dynamics of retrotransposition-driven genomic expansions in Oryza australiensis, a wild relative of rice. Genome Res. 2006;16:1262–9. Bennetzen JL, Kellogg EA. Do Plants Have a One-Way Ticket to Genomic Obesity? Plant Cell [Internet]. 1997;9:1509–14. Available from: http://www.plantcell.org/content/9/9/1509.short Munasinghe M, Read A, Stitzer M, Song B, Menard C, Ma K, et al. Combined analysis of transposable elements and structural variation in maize genomes reveals genome contraction outpaces expansion. PLoS Genet. 2023;19:e1011086. Devos KM, Brown JKM, Bennetzen JL. Genome size reduction through illegitimate recombination counteracts genome expansion in Arabidopsis. Genome Res. 2002;12:1075–9. Tian Z, Rizzon C, Du J, Zhu L, Bennetzen JL, Jackson SA, et al. Do genetic recombination and gene density shape the pattern of DNA elimination in rice long terminal repeat retrotransposons? Genome Res. 2009;19:2221–30. Shirasu K, Schulman A, Lahaye T, Schulze-Lefert P. A contiguous 66-kb barley DNA sequence provides evidence for reversible genome expansion. Genome Res. 2000;10:908–15. Vitte C, Panaud O. Formation of solo-LTRs through unequal homologous recombination counterbalances amplifications of LTR retrotransposons in rice Oryza sativa L. Mol Biol Evol. 2003;20(4):528 – 40. Mol Biol Evol. 2003;20:528–40. Viot CR, Wendel JF. Evolution of the Cotton Genus, Gossypium, and Its Domestication in the Americas. CRC Crit Rev Plant Sci [Internet]. 2023;42:1–33. Available from: https://doi.org/10.1080/07352689.2022.2156061 Daojun Y, Grover CE, Hu G, Pan M, Miller ER, Conover JL, et al. Parallel and Intertwining Threads of Domestication in Allopolyploid Cotton. Adv Sci. 2021;8:2003634. Wang M, Li J, Wang P, Liu F, Liu Z, Zhao G, et al. Comparative Genome Analyses Highlight Transposon-Mediated Genome Expansion and the Evolutionary Architecture of 3D Genomic Folding in Cotton. Mol Biol Evol. 2021;38:3621–36. He X, Qi Z, Liu Z, Chang X, Zhang X, Li J et al. Pangenome analysis reveals transposon-driven genome evolution in cotton. BMC Biol [Internet]. 2024;22:92. Available from: https://doi.org/10.1186/s12915-024-01893-2 Chang X, He X, Li J, Liu Z, Pi R, Luo X, et al. High-quality Gossypium hirsutum and Gossypium barbadense genome assemblies reveal the landscape and evolution of centromeres. Plant Commun. 2024;5:100722. Li J, Liu Z, You C, Qi Z, You J, Grover CE et al. Convergence and divergence of diploid and tetraploid cotton genomes. Nat Genet. 2024;10.1038/s41588-024-01964–8. Naish M, Henderson IR. The structure, function, and evolution of plant centromeres. Genome Res. 2024;34:161–78. Yi X, Liang Y, Huerta-Sanchez E, Jin1 X, Cuo Z, Ping X et al. Sequencing of Fifty Human Exomes Reveals Adaptation to High Altitude. Science (1979). 2010;329:75–78. Hu Y, Chen J, Fang L, Zhang Z, Ma W, Niu Y et al. Gossypium barbadense and Gossypium hirsutum genomes provide insights into the origin and evolution of allotetraploid cotton. Nature Genetics. 2019 51:4 [Internet]. 2019 [cited 2025 Feb 21];51:739–48. Available from: https://www.nature.com/articles/s41588-019-0371-5 Bousios A, Kakutani T, Henderson IR. Centrophilic Retrotransposons of Plant Genomes. Annu Rev Plant Biol [Internet]. 2025 [cited 2025 Mar 6]; Available from: https://pubmed.ncbi.nlm.nih.gov/39952673/ Schneider KL, Xie Z, Wolfgruber TK, Presting GG. Inbreeding drives maize centromere evolution. Proc Natl Acad Sci U S A [Internet]. 2016 [cited 2025 Feb 14];113:E987–96. Available from: https://pubmed.ncbi.nlm.nih.gov/26858403/ Wang Z, Wang W, Xie X, Wang Y, Yang Z, Peng H et al. Dispersed emergence and protracted domestication of polyploid wheat uncovered by mosaic ancestral haploblock inference. Nat Commun [Internet]. 2022 [cited 2025 Feb 14];13:3891. Available from: https://pmc.ncbi.nlm.nih.gov/articles/PMC9259585/ Meyer RS, Purugganan MD. Evolution of crop species: Genetics of domestication and diversification. Nat Rev Genet [Internet]. 2013;14:840–52. Available from: http://www.ncbi.nlm.nih.gov/pubmed/24240513 Grover CE, Yoo MJ, Lin M, Murphy MD, Harker DB, Byers RL et al. Genetic Analysis of the Transition from Wild to Domesticated Cotton (Gossypium hirsutum L.). G3: Genes|Genomes|Genetics [Internet]. 2019 [cited 2025 Feb 21];10:731. Available from: https://pmc.ncbi.nlm.nih.gov/articles/PMC7003101/ Teale WD, Paponov IA, Palme K. Auxin in action: signalling, transport and the control of plant growth and development. Nature Reviews Molecular Cell Biology 2006 7:11 [Internet]. 2006 [cited 2025 Feb 21];7:847–59. Available from: https://www.nature.com/articles/nrm2020 Jareczek JJ, Grover CE, Wendel JF. Cotton fiber as a model for understanding shifts in cell development under domestication. Front Plant Sci [Internet]. 2023;14. Available from: https://www.frontiersin.org/journals/plant-science/articles/ 10.3389/fpls.2023.1146802 Chu Q, Fu X, Zhao J, Li Y, Liu L, Zhang L et al. Simultaneous improvement of fiber yield and quality in upland cotton (Gossypium hirsutum L.) by integration of auxin transport and synthesis. Molecular Breeding [Internet]. 2024 [cited 2025 Feb 21];44:1–17. Available from: https://link.springer.com/article/10.1007/s11032-024-01500-w Ma J, Pei W, Ma Q, Geng Y, Liu G, Liu J et al. QTL analysis and candidate gene identification for plant height in cotton based on an interspecific backcross inbred line population of Gossypium hirsutum × Gossypium barbadense. Theoretical and Applied Genetics [Internet]. 2019 [cited 2025 Feb 21];132:2663–76. Available from: https://link.springer.com/article/ 10.1007/s00122-019-03380-7 Shi H, Qanmber G, Yang Z, Guo Y, Ma S, Shu S, et al. An AP2/ERF transcription factor GhERF109 negatively regulates plant growth and development in cotton. Plant Sci. 2025;352:112365. Stortenbeker N, Bemer M. The SAUR gene family: the plant’s toolbox for adaptation of growth and development. J Exp Bot [Internet]. 2019 [cited 2025 Feb 24];70:17–27. Available from: https://pubmed.ncbi.nlm.nih.gov/30239806/ Gonzalez-Ibeas D, Ibanez V, Perez-Roman E, Borredá C, Terol J, Talon M. Shaping the biology of citrus: II. Genomic determinants of domestication. Plant Genome [Internet]. 2021 [cited 2025 Feb 24];14. Available from: https://pubmed.ncbi.nlm.nih.gov/34464512/ Zhao Z, Chen T, Yue J, Pu N, Liu J, Luo L et al. Small Auxin Up RNA 56 (SAUR56) regulates heading date in rice. Mol Breed [Internet]. 2023 [cited 2025 Feb 24];43. Available from: https://pubmed.ncbi.nlm.nih.gov/37521314/ Guo N, Wang S, Wang T, Duan M, Zong M, Miao L, et al. A graph-based pan-genome of Brassica oleracea provides new insights into its domestication and morphotype diversification. Plant Commun. 2024;5:100791. Ou S, Su W, Liao Y, Chougule K, Agda JRA, Hellinga AJ, et al. Benchmarking transposable element annotation methods for creation of a streamlined, comprehensive pipeline. Genome Biol. 2019;20:275. Zhang R-G, Li G-Y, Wang X-L, Dainat J, Wang Z-X, Ou S et al. TEsorter: An accurate and fast method to classify LTR-retrotransposons in plant genomes. Hortic Res [Internet]. 2022;9:uhac017. Available from: https://doi.org/10.1093/hr/uhac017 Fu L, Niu B, Zhu Z, Wu S, Li W. CD-HIT: Accelerated for clustering the next-generation sequencing data. Bioinformatics. 2012;28:3150–2. Ou S, Jiang N, LTR_retriever:. A Highly Accurate and Sensitive Program for Identification of Long Terminal Repeat Retrotransposons. Plant Physiol [Internet]. 2018 [cited 2025 Jan 7];176:1410–22. Available from: https://pubmed.ncbi.nlm.nih.gov/29233850/ Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403–10. Li H. Minimap2: Pairwise alignment for nucleotide sequences. Bioinformatics. 2018;34:3094–100. Heller D, Vingron M. SVIM-asm: structural variant detection from haploid and diploid genome assemblies. Bioinformatics [Internet]. 2021;36:5519–21. Available from: https://doi.org/10.1093/bioinformatics/btaa1034 Danecek P, Bonfield JK, Liddle J, Marshall J, Ohan V, Pollard MO et al. Twelve years of SAMtools and BCFtools. Gigascience [Internet]. 2021;10:giab008. Available from: https://doi.org/10.1093/gigascience/giab008 English AC, Menon VK, Gibbs RA, Metcalf GA, Sedlazeck FJ. Truvari: refined structural variant comparison preserves allelic diversity. Genome Biol [Internet]. 2022;23:271. Available from: https://doi.org/10.1186/s13059-022-02840-6 Hickey G, Heller D, Monlong J, Sibbesen JA, Sirén J, Eizenga J et al. Genotyping structural variants in pangenome graphs using the vg toolkit. Genome Biol [Internet]. 2020;21:35. Available from: https://doi.org/10.1186/s13059-020-1941-7 Granato ISC, Galli G, de Oliveira Couto EG, e Souza MB, Mendonça LF, Fritsche-Neto R. snpReady: a tool to assist breeders in genomic analysis. Molecular Breeding [Internet]. 2018 [cited 2025 Jan 7];38:1–7. Available from: https://link.springer.com/article/ 10.1007/s11032-018-0844-8 Klopfenstein DV, Zhang L, Pedersen BS, Ramírez F, Vesztrocy AW, Naldi A et al. GOATOOLS: A Python library for Gene Ontology analyses. Scientific Reports 2018 8:1 [Internet]. 2018 [cited 2025 Jan 7];8:1–17. Available from: https://www.nature.com/articles/s41598-018-28948-z Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29:644–52. Patro R, Duggal G, Love MI, Irizarry RA, Kingsford C. Salmon provides fast and bias-aware quantification of transcript expression. Nat Methods [Internet]. 2017 [cited 2025 Feb 27];14:417–9. Available from: https://pubmed.ncbi.nlm.nih.gov/28263959/ Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550. Additional Declarations No competing interests reported. Supplementary Files supfigs.pdf STables.xlsx Cite Share Download PDF Status: Published Journal Publication published 27 Oct, 2025 Read the published version in Genome Biology → Version 1 posted Editorial decision: Revision requested 21 May, 2025 Reviews received at journal 05 May, 2025 Reviews received at journal 02 May, 2025 Reviews received at journal 26 Apr, 2025 Reviewers agreed at journal 18 Apr, 2025 Reviewers agreed at journal 15 Apr, 2025 Reviews received at journal 11 Apr, 2025 Reviewers agreed at journal 01 Apr, 2025 Reviewers agreed at journal 01 Apr, 2025 Reviewers agreed at journal 31 Mar, 2025 Reviewers agreed at journal 31 Mar, 2025 Reviewers invited by journal 31 Mar, 2025 Editor assigned by journal 17 Mar, 2025 Submission checks completed at journal 10 Mar, 2025 First submitted to journal 06 Mar, 2025 You are reading this latest preprint version Research Square lets you share your work early, gain feedback from the community, and start making changes to your manuscript prior to peer review in a journal. As a division of Research Square Company, we’re committed to making research communication faster, fairer, and more useful. We do this by developing innovative software and high quality services for the global research community. Our growing team is made up of researchers and industry professionals working together to solve the most critical problems facing scientific publishing. Also discoverable on Platform About Our Team In Review Editorial Policies Advisory Board Help Center Resources Author Services Accessibility API Access RSS feed Manage Cookie Preferences © Research Square 2026 | ISSN 2693-5015 (online) Privacy Policy Terms of Service Do Not Sell My Personal Information {"props":{"pageProps":{"initialData":{"identity":"rs-6172192","acceptedTermsAndConditions":true,"allowDirectSubmit":false,"archivedVersions":[],"articleType":"Research Article","associatedPublications":[],"authors":[{"id":429872259,"identity":"fffce823-509b-434e-b833-593f16a3505a","order_by":0,"name":"Lucía Campos-Dominguez","email":"","orcid":"","institution":"CRAG (CSIC- IRTA-UAB-UB)","correspondingAuthor":false,"prefix":"","firstName":"Lucía","middleName":"","lastName":"Campos-Dominguez","suffix":""},{"id":429872262,"identity":"e9e74556-7c82-49de-9843-efd2d6919cc7","order_by":1,"name":"Raúl Castanera","email":"","orcid":"","institution":"CRAG (CSIC- IRTA-UAB-UB)","correspondingAuthor":false,"prefix":"","firstName":"Raúl","middleName":"","lastName":"Castanera","suffix":""},{"id":429872264,"identity":"2328f281-4550-466a-81dd-6afe9892c239","order_by":2,"name":"Corrinne E. Grover","email":"","orcid":"","institution":"Iowa State University","correspondingAuthor":false,"prefix":"","firstName":"Corrinne","middleName":"E.","lastName":"Grover","suffix":""},{"id":429872266,"identity":"daa4ed08-fc47-4045-b5db-b2b9232fcf4a","order_by":3,"name":"Jonathan F. Wendel","email":"","orcid":"","institution":"Iowa State University","correspondingAuthor":false,"prefix":"","firstName":"Jonathan","middleName":"F.","lastName":"Wendel","suffix":""},{"id":429872268,"identity":"3ad9f627-741e-4430-bad4-ed685493286f","order_by":4,"name":"Josep M. Casacuberta","email":"data:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAZAAAAAyAQMAAABI0h/eAAAABlBMVEX///8AAABVwtN+AAAACXBIWXMAAA7EAAAOxAGVKw4bAAABR0lEQVRIiWNgGAWjYBACxgYkDjOYZG/AKolbiwQDzwH8WlAARItEAgFV7e0PHxcw2MnLtx9+Jl1Qca+Of+bjZx9+/LLJ4xc73PjhB4ONPbrDes4YG89gSDbccCbNTHrGmWIJidtpxjN7+9KKJWcnNkv2MKQlojmPcUYOmzQPwwHGDQxABm9bggTD7RxmBt6ew4kbbie2MfAwHEZ3J+OM9Oe/gVrs5/e/AWr5lyAhf/MMM+Pfnv9gLYx/GP5jOGxGghkzUEtiww2QLQ0JEgY3eJiZeX4cAGsBSaGHG8gv0jMMkpM33HhmbD3jWILkxjNpxsyyDcmJM4F+kZYxSEb3iyEwxD4XVNjZzu9Pfni7oCaBX+744ceMb/7YJfZLpz/8+KbCDt1hhg2g6DBAD3vGNhgLQ4pBngGWTlDBHyxio2AUjIJRMFIBABXxcdV7hnRGAAAAAElFTkSuQmCC","orcid":"","institution":"CRAG (CSIC- IRTA-UAB-UB)","correspondingAuthor":true,"prefix":"","firstName":"Josep","middleName":"M.","lastName":"Casacuberta","suffix":""}],"badges":[],"createdAt":"2025-03-06 16:08:19","currentVersionCode":1,"declarations":"","doi":"10.21203/rs.3.rs-6172192/v1","doiUrl":"https://doi.org/10.21203/rs.3.rs-6172192/v1","draftVersion":[],"editorialEvents":[{"content":"https://doi.org/10.1186/s13059-025-03837-7","type":"published","date":"2025-10-27T15:58:08+00:00"}],"editorialNote":"","failedWorkflow":false,"files":[{"id":78937081,"identity":"74ffdc94-200c-4f3c-89d1-b86cf48667fa","added_by":"auto","created_at":"2025-03-21 05:34:34","extension":"png","order_by":1,"title":"Figure 1","display":"","copyAsset":false,"role":"figure","size":49140,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eLTR-RT content in diploid and tetraploid cottons by subgenome.\u003c/strong\u003e \u003cstrong\u003ea) \u003c/strong\u003eNumber of intact LTR-RT per megabase pairs per subgenome in each of the cotton species and varieties analyzed (left panel), and violin plots representing the insertion times intact LTR-RTs per subgenome in each of the cotton species and varieties analyzed (right panel). \u003cstrong\u003eb)\u003c/strong\u003e Lineage classification of the intact LTR-RT elements of each species and subgenome. The values of the tetraploid \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e species are represented by the average values of all varieties of the corresponding species analysed (left panel). Total assembled genome % corresponding to each LTR-RT lineage in the diploid and tetraploid cotton species analysed. The tetraploid \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e genomes are represented by the varieties TM1 and Pima90, respectively (right panel).\u003c/p\u003e","description":"","filename":"Fig1.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/d95c59b3ec29fb4c6106b762.png"},{"id":78937083,"identity":"4d6f6ab9-0a7b-4675-9c0f-a9e10e9a08fc","added_by":"auto","created_at":"2025-03-21 05:34:34","extension":"png","order_by":2,"title":"Figure 2","display":"","copyAsset":false,"role":"figure","size":65500,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eIntact LTR-RT and solo LTR \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eCopia \u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003eTork and \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eGypsy\u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003e Athila content in diploid and tetraploid cottons by subgenome\u003c/strong\u003e. \u0026nbsp;Number of intact LTR-RT per megabase pairs per subgenome in each of the cotton species and varieties analyzed (top panels), violin plots representing the insertion times these intact LTR-RTs per subgenome in each of the cotton species and varieties analyzed (middle panels), and number of solo LTRs per megabase pairs per subgenome in each of the cotton species and varieties analyzed (bottom panels).\u003c/p\u003e","description":"","filename":"Fig2.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/1f236c7e3f2689cff48c3a81.png"},{"id":78937718,"identity":"055ac6b7-35bf-4976-87b8-86e08b607dd2","added_by":"auto","created_at":"2025-03-21 05:42:34","extension":"png","order_by":3,"title":"Figure 3","display":"","copyAsset":false,"role":"figure","size":69287,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eIntact LTR-RT and solo LTR \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eGypsy \u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003eCRM and \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eGypsy\u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003e Athila content in diploid and tetraploid cottons by subgenome\u003c/strong\u003e. \u0026nbsp;Number of intact LTR-RT per megabase pairs per subgenome in each of the cotton species and varieties analyzed (top panels), violin plots representing the insertion times these intact LTR-RTs per subgenome in each of the cotton species and varieties analyzed (middle panels), and number of solo LTRs per megabase pairs per subgenome in each of the cotton species and varieties analyzed (bottom panels).\u003c/p\u003e","description":"","filename":"Fig3.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/49219c7c845c0b077b0fb4a6.png"},{"id":78937080,"identity":"2e643bac-c5e5-438a-9a5b-25b7b0077edc","added_by":"auto","created_at":"2025-03-21 05:34:34","extension":"png","order_by":4,"title":"Figure 4","display":"","copyAsset":false,"role":"figure","size":45785,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eChromosomal distribution of CRM and Tekay elements.\u003c/strong\u003e \u003cstrong\u003ea) \u003c/strong\u003eDensity of CRM and Tekay intact LTR-RTs and solo-LTRs along chromosome D07 of \u003cem\u003eG. hirsutum\u003c/em\u003e. \u003cstrong\u003eb) \u003c/strong\u003eDensity of CRM and Tekay intact LTR-RTs and solo-LTRs along chromosomes A10 and D07 in the diploid and tetraploid cottons here analysed. A green asterisk indicates where the scale of the x axis of solo-LTRs represented has been modified (2x) in order allow visualizing their distribution along chromosomes.\u003c/p\u003e","description":"","filename":"Fig4.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/2f9c1f458ace9a08d6326d72.png"},{"id":78937082,"identity":"8339f237-f87e-46ff-a1e8-9d1390fb8515","added_by":"auto","created_at":"2025-03-21 05:34:34","extension":"png","order_by":5,"title":"Figure 5","display":"","copyAsset":false,"role":"figure","size":31783,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eIntraspecific variability of LTR-RT and solo-LTR insertions.\u003c/strong\u003e \u0026nbsp;\u003cstrong\u003ea)\u003c/strong\u003e Number of TIPs corresponding to intact LTR-RT and solo-LTR inserions in the \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e pangenomes per LTR-RT lineage. \u003cstrong\u003eb)\u003c/strong\u003e Percentage of Tork and Ivana TIPs in different genome regions as compared with a distribution at random (percentage of the genome occupied by each region). Asterisks represent significant enrichment/depletion of insertions in each region compared to random insertions (p-value \u0026lt; 0.05).\u003c/p\u003e","description":"","filename":"Fig5.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/4e53f3dccd3f9142c71972aa.png"},{"id":78937086,"identity":"38c3ac41-3532-455b-8f3e-9bf64c64eb53","added_by":"auto","created_at":"2025-03-21 05:34:35","extension":"png","order_by":6,"title":"Figure 6","display":"","copyAsset":false,"role":"figure","size":73299,"visible":true,"origin":"","legend":"\u003cp\u003e\u003cstrong\u003eChanges in TIP frequencies accompanying \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eG. hirsutum\u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003e and \u003c/strong\u003e\u003cem\u003e\u003cstrong\u003eG. barbadense\u003c/strong\u003e\u003c/em\u003e\u003cstrong\u003e domestication and breeding.\u003c/strong\u003e \u003cstrong\u003ea)\u003c/strong\u003e PCA based on TIPs in wild, landrace and domesticated varieties of \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e. \u003cstrong\u003eb)\u003c/strong\u003e Distribution of Population Branch Statistic values for TIPs (deletions and insertions with respect to the reference genome) for the comparison of cultivated against landrace and wild accessions in \u003cem\u003eG. hirsutum.\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e. \u0026nbsp;\u003cstrong\u003ec)\u003c/strong\u003ePercentage of auxin-response GO term in the group of genes linked to a TIP with high PBS (5% top PBS) compared with the percentage of this GO term in the whole \u003cem\u003eG. hirsutum\u003c/em\u003e proteome. \u003cstrong\u003ed)\u003c/strong\u003e Fst-based tree representing the distance between the three populations using all shared TIPs (n= 5,535), a deletion TIP 565 bp upstream of a SAUR-like gene (Gh_D13G264300, auxin-responsive gene) potentially related with domestication (D13:63,658,839) and a deletion TIP present in the intron of a SAUR-like gene (Gh_A09G167200, auxin-responsive gene) potentially related with cotton diversification.\u003c/p\u003e","description":"","filename":"Fig6.png","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/087b21116e8e9f60d09cad43.png"},{"id":95040458,"identity":"dbaaa286-34d8-48a2-83f7-c1c64813e68b","added_by":"auto","created_at":"2025-11-03 16:09:00","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"manuscript-pdf","size":1336691,"visible":true,"origin":"","legend":"","description":"","filename":"manuscript.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/7524a4a8-81a4-494e-9223-f2561f5a9efe.pdf"},{"id":78937099,"identity":"40dded20-eedb-4cac-9450-96357b234a52","added_by":"auto","created_at":"2025-03-21 05:34:35","extension":"pdf","order_by":0,"title":"","display":"","copyAsset":false,"role":"supplement","size":11503481,"visible":true,"origin":"","legend":"","description":"","filename":"supfigs.pdf","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/aaeedc80aeb33efeae6f3921.pdf"},{"id":78937721,"identity":"87f2edfe-c9dd-4132-8a7a-e0e0cf9f882d","added_by":"auto","created_at":"2025-03-21 05:42:35","extension":"xlsx","order_by":1,"title":"","display":"","copyAsset":false,"role":"supplement","size":111763,"visible":true,"origin":"","legend":"","description":"","filename":"STables.xlsx","url":"https://assets-eu.researchsquare.com/files/rs-6172192/v1/1b39842b2ba4320cb5247133.xlsx"}],"financialInterests":"No competing interests reported.","formattedTitle":"Differential LTR-retrotransposon dynamics across polyploidization, speciation, domestication and improvement of cotton (Gossypium)","fulltext":[{"header":"Background","content":"\u003cp\u003eTransposable elements (TEs) are major components of plant genomes. Their repetitive nature and their ability to amplify and move within the genome makes them a major source of genetic variability and one of the major drivers of plant genome evolution [\u003cspan citationid=\"CR1\" class=\"CitationRef\"\u003e1\u003c/span\u003e, \u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e]. TEs are an important source of mutations underlying crop domestication and improvement [\u003cspan citationid=\"CR3\" class=\"CitationRef\"\u003e3\u003c/span\u003e]. However, as potential mutagens TEs are tightly controlled, mainly through epigenetic mechanisms [\u003cspan citationid=\"CR4\" class=\"CitationRef\"\u003e4\u003c/span\u003e]. For this reason, the dynamics of TEs within genomes has been often seen as an arms-race between an invasive parasite and a defensive host, with genomes tightly controlling TE activity and TEs evolving to evade genomic suppression [\u003cspan additionalcitationids=\"CR6\" citationid=\"CR5\" class=\"CitationRef\"\u003e5\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR7\" class=\"CitationRef\"\u003e7\u003c/span\u003e]. TE activity is known to increase in situations where DNA methylation and epigenetic silencing decrease, such as stress situations [\u003cspan citationid=\"CR8\" class=\"CitationRef\"\u003e8\u003c/span\u003e], including the genomic shock that results from the combination of two different genomes in polyploid species [\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e, \u003cspan citationid=\"CR10\" class=\"CitationRef\"\u003e10\u003c/span\u003e]. This activation can lead to a burst of transposition and a subsequent repression when silencing is re-established [\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e]. Yet accumulated data during the last few years in plant genomes suggests that the dynamics of different TE types, lineages and families that coexist in a genome are highly diverse. These lineages and families may thus be viewed as different players in an evolving ecosystem [\u003cspan additionalcitationids=\"CR13\" citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eLong terminal repeat retrotransposons (LTR-RTs) are the most prevalent type of TEs in plants, as exemplified by maize, where LTR-RTs account for 90% of the TEs and more than 80% of the genome space [\u003cspan citationid=\"CR14\" class=\"CitationRef\"\u003e14\u003c/span\u003e]. LTR-RTs are divided into two major superfamilies, \u003cem\u003eCopia\u003c/em\u003e and \u003cem\u003eGypsy\u003c/em\u003e, which often show important differences in their regulation and distribution within genomes [\u003cspan citationid=\"CR15\" class=\"CitationRef\"\u003e15\u003c/span\u003e]. LTR-RTs transpose through a replicative mechanism, i.e., via an RNA copy reverse transcribed into cDNA prior to integration at a new location. Therefore, their activity results in an increase in copy number, which endows them with the potential to invade, and thereby inflate, the genome. As an example, the genome of the wild rice \u003cem\u003eO. australiensis\u003c/em\u003e doubled in size in just three million years by the amplification of only three families of LTR-RTs [\u003cspan citationid=\"CR16\" class=\"CitationRef\"\u003e16\u003c/span\u003e]. Counterbalancing this genome size expansion via TE proliferation, LTR-RTs can also be eliminated from genomes by different mechanisms [\u003cspan citationid=\"CR17\" class=\"CitationRef\"\u003e17\u003c/span\u003e, \u003cspan citationid=\"CR18\" class=\"CitationRef\"\u003e18\u003c/span\u003e]. The main mechanism for this deletional process is illegitimate recombination, often at the LTRs, giving rise to so-called solo-LTRs resulting in truncated LTR-copies [\u003cspan additionalcitationids=\"CR20 CR21\" citationid=\"CR19\" class=\"CitationRef\"\u003e19\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR22\" class=\"CitationRef\"\u003e22\u003c/span\u003e]. Although the importance of recombination-mediated LTR-RT elimination is well recognized, its analysis is often neglected when studying LTR-RT dynamics.\u003c/p\u003e \u003cp\u003eThe cotton genus (\u003cem\u003eGossypium\u003c/em\u003e) serves as an excellent evolutionary model for polyploidization, domestication and crop improvement. Allopolyploid cotton species arose 1\u0026ndash;2 mya from the interspecific hybridisation of an \u0026ldquo;A-genome\u0026rdquo; and a \u0026ldquo;D-genome\u0026rdquo; species, which diverged from each other 5\u0026ndash;7 mya [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e]. Because of this long divergence of the two genome donors, A and D, their resulting subgenomes can be easily distinguished in the tetraploid genomes. The ancestor providing the A subgenome is closely related to the modern species \u003cem\u003eG\u003c/em\u003e. \u003cem\u003eherbaceum\u003c/em\u003e (A1) and \u003cem\u003eG. arboreum\u003c/em\u003e (A2), whereas the ancestor of the D subgenome is closest, among extant D-genome species, to \u003cem\u003eG. raimondii\u003c/em\u003e (D5) [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e]. Two allotetraploid species (\u003cem\u003eG. hirsutum\u003c/em\u003e, AD1, and \u003cem\u003eG. barbadense\u003c/em\u003e, AD2) were independently domesticated from different wild progenitors in Central and South America, respectively, over the last 8000 years [\u003cspan citationid=\"CR23\" class=\"CitationRef\"\u003e23\u003c/span\u003e, \u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eTE studies on cotton species have previously highlighted their significant role on genome size expansion, being the major contributors to genome size differences between diploid cotton species [\u003cspan citationid=\"CR25\" class=\"CitationRef\"\u003e25\u003c/span\u003e]. Recent analyses have also suggested that TEs may have amplified accompanying the polyploidization event, particularly for the \u003cem\u003eGypsy\u003c/em\u003e Tekay and CRM lineages of LTR-RTs [\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e], which may have reshaped the centromeres in tetraploids [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eHere we study the LTR-RT dynamics (both LTR-RT insertion and removal through solo-LTR formation) accompanying diploid divergence, polyploidization, and speciation at the allopolyploid level using genome assemblies of eight \u003cem\u003eG. hirsutum\u003c/em\u003e and ten \u003cem\u003eG. barbadense\u003c/em\u003e accessions and representatives of their parental diploid species. Using pangenomes generated for the two domesticated tetraploid cotton species (\u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e), we report on recent LTR-RT dynamics accompanying initial crop domestication and subsequent improvement. Using this temporally stratified and explicit approach, we show that each LTR-RT lineage has experienced its own dynamics throughout the different events characterizing the recent evolution of cotton species and varieties. We also use the pangenomes to analyze relatively recent or ongoing transposition events within the two allopolyploid species, documenting the frequency of the intraspecies transposon insertion polymorphisms (TIPs) in approximately 500 sequenced genomes representing the wild-to-domesticated continuum within \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense.\u003c/em\u003e Our results show that LTR-RTs may have been instrumental in shaping cotton genomes and generating the genome variability that was subject to selection during domestication and crop improvement.\u003c/p\u003e"},{"header":"Results","content":"\u003cdiv id=\"Sec3\" class=\"Section2\"\u003e \u003ch2\u003eLTR-RT lineages have differentially amplified in tetraploid cotton species\u003c/h2\u003e \u003cp\u003eAs a first step in exploring the dynamics of LTR-RTs in the recent evolution of cotton, we annotated intact LTR-RT copies in the highly contiguous assemblies of eight accessions of \u003cem\u003eG. hirsutum\u003c/em\u003e and ten accessions of \u003cem\u003eG. barbadense\u003c/em\u003e, as well as in the genomes of two diploid species that represent the best models of the ancestral A and D subgenomes of cotton (\u003cem\u003eG. herbaceum\u003c/em\u003e and \u003cem\u003eG. raimondii\u003c/em\u003e, respectively). Our analysis focused only on the intact LTR-RT copies (identified by stringent approaches; see methods) to concentrate on the youngest elements, i.e., those that could have arisen from transposition accompanying polyploidization, speciation, and domestication of cotton, events that occurred within the last 1\u0026ndash;2\u0026nbsp;million years.\u003c/p\u003e \u003cp\u003eResults in Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea (left) show that the two polyploid species have, in general, a higher number of intact LTR-RTs per megabase than do the diploids, especially in the D subgenome, where the density of intact LTR-RT elements is roughly doubled (from 3.5 elements per Mbp in the D-genome diploid parental \u003cem\u003eG. raimondii\u003c/em\u003e to ~\u0026thinsp;8 and ~\u0026thinsp;7 elements per Mbp in \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e, respectively). For the A subgenome, there is an increased density with respect to that of \u003cem\u003eG. herbaceum\u003c/em\u003e but not with respect to \u003cem\u003eG. arboreum\u003c/em\u003e. However, it has been shown that \u003cem\u003eG. arboreum\u003c/em\u003e experienced a recent increase in LTR-RT copy number, and therefore its present LTR-RT density does not reflect its density when polyploidization occurred [\u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e]. Therefore, we only use \u003cem\u003eG. herbaceum\u003c/em\u003e as a comparator for the A subgenomes henceforth. Between the two tetraploid species, \u003cem\u003eG. hirsutum\u003c/em\u003e varieties generally have a slightly higher number of intact LTR-RTs than the \u003cem\u003eG. barbadense\u003c/em\u003e accessions, and there are no major differences among accessions within each of the two species (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea, left panel).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eInsertion times were calculated for each element by comparing the two LTR sequences for each intact LTR-RT. Figure\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea (right) shows that LTR-RT insertions are younger, on average, in the polyploids than in the diploids (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003ea, right panel), suggesting more recent activation of LTR-RTs in the polyploids. Analysis of the LTR-RT composition of the genome of the two tetraploid species relative to their two diploid parentals reveals that the increase in intact copies, and therefore in the percentage of the genome covered by LTR-RTs, is mainly due to few LTR-RT lineages (e.g., CRM and Tork, see Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003eAnalysis of the number of intact elements from each LTR-RT lineage (Figs.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e, \u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e and Supplementary Fig.\u0026nbsp;1) indicates that there is diversity in dynamics even among lineages belonging to the same superfamily of LTR-RTs (i.e. \u003cem\u003eCopia\u003c/em\u003e and \u003cem\u003eGypsy\u003c/em\u003e), whereby some lineages exhibit abundant increases in copy number in the polyploids whereas others remain relatively stable during recent cotton evolution. For example, copy number for the \u003cem\u003eCopia\u003c/em\u003e Tork lineage is three times higher per megabase in the two tetraploids than in the two diploid progenitors (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e top left). Insertion time analysis of Tork elements shows that although most elements are relatively old in the diploids, they are notably younger in the polyploids (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e middle left). This confirms that Tork has been particularly active in the tetraploid genomes and may still be actively transposing. In contrast, the number of copies per megabase of other lineages, such as the relatively abundant \u003cem\u003eGypsy\u003c/em\u003e Athila (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e top right), is similar when comparing each subgenome of the polyploid with its corresponding diploid genome, suggesting no or little active retrotransposition since polyploid formation. An analysis of the insertion time of Athila intact elements shows a relatively old insertion time, particularly for elements residing in the A subgenomes, and an almost identical distribution of insertion times between diploid and tetraploid genomes, supporting very low activity after polyploidization (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e middle right). However, Athila elements have probably been very active in the past, as Athila-related sequences account for a sizeable fraction (~\u0026thinsp;13% and 5% in subgenomes A and D, respectively) of the genome of both tetraploid species and their diploid parentals (~\u0026thinsp;10% in \u003cem\u003eG. herbaceum\u003c/em\u003e and 4,6% in \u003cem\u003eG. raimondii\u003c/em\u003e) (Fig.\u0026nbsp;\u003cspan refid=\"Fig1\" class=\"InternalRef\"\u003e1\u003c/span\u003eb).\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eA general inference from the analysis of copy number for all LTR-RT lineages is that the evolutionary behavior of the different TE families varies considerably among lineages, subgenomes, species and, with respect to timing of proliferation events (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e, \u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e and Supplementary Fig.\u0026nbsp;1). Some LTR-RT lineages show an increase in retrotransposition in both subgenomes after polyploidisation (the \u003cem\u003eCopias\u003c/em\u003e Tork and Ivana, and the \u003cem\u003eGypsy\u003c/em\u003e CRM), whereas others mainly increased in one subgenome (the \u003cem\u003eGypsy\u003c/em\u003e Tekay). In contrast, a stable number of LTR-RT insertions was observed among all species at both the diploid and tetraploid level for some retrotransposon lineages, such as the \u003cem\u003eGypsy\u003c/em\u003e Athila, Ogre, Reina, Galadriel families and the \u003cem\u003eCopia\u003c/em\u003e Ale family. Finally, some TE families are present at similarly low copy numbers in all genomes (the \u003cem\u003eCopias\u003c/em\u003e Sire, Bianca, Ikeros and TAR). In general, lineages present at a relatively high copy number exhibited similar insertion numbers when comparing different varieties of the same species, although this is not always the case, as evidenced by the Ivana lineage in the subgenome D of \u003cem\u003eG. hirsutum\u003c/em\u003e (Supplementary Fig.\u0026nbsp;1). Similarly, the copy number of each lineage is generally similar between the two tetraploid species, although differences for specific lineages do exist. Notably among these exceptions are the CRM and Tekay \u003cem\u003eGypsy\u003c/em\u003e lineages (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e), which are further discussed below.\u003c/p\u003e \u003c/div\u003e\n\u003ch3\u003eThe rate of LTR-RT elimination varies among lineages, species and subgenomes\u003c/h3\u003e\n\u003cp\u003eWe further explored the dynamics of LTR-RT elements accompanying recent cotton evolution by quantifying LTR-RT elimination though intra element LTR-recombination. We annotated potential solo-LTRs by searching for LTRs corresponding to the lineages and families of intact LTR-RTs that did not have surrounding sequences that could correspond to LTR-RT internal sequences (see methods). These results show that there are important differences among the LTR-RT lineages (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e, \u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e and supplementary Fig.\u0026nbsp;1 bottom panels). For example, the relative number of Tork solo-LTRs is similar between tetraploid and diploid genomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e bottom left), despite the increase in intact Tork insertions in polyploid cotton relative to their parental diploids (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e top left). Combined with the younger age of intact Tork elements in the tetraploids (relative to the diploids; Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e middle left), these observations suggest both that Tork elements have been more active in the tetraploids than in the diploids (after polyploidization) and that only a minor fraction of the new insertions has been removed by recombination. In contrast, both diploids and tetraploids contain two to ten times more Athila solo-LTRs than intact elements (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e right), suggesting a high rate of elimination of these elements. As noted above, Athila has similarly old element insertion times between diploids and polyploids genomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig2\" class=\"InternalRef\"\u003e2\u003c/span\u003e middle right), suggesting low or no recent activity of this family; however, the greater number of solo-LTRs in \u003cem\u003eG. barbadense\u003c/em\u003e relative to \u003cem\u003eG. hirsutum\u003c/em\u003e may indicate recent activity and rapid turnover of Athila insertions in \u003cem\u003eG. barbadense\u003c/em\u003e.\u003c/p\u003e \u003cp\u003eA similar diversity of solo-LTR dynamics characterizes the remaining LTR-RT lineages. For example, Ogre solo-LTRs are also more abundant than intact elements (supplementary Fig.\u0026nbsp;1); however, unlike Athila, intact elements tend to be younger in the tetraploids, particularly for those in the D subgenome of \u003cem\u003eG. hirsutum\u003c/em\u003e, where the number of solo-LTRs is particularly high. This observation may suggest that the removal of intact elements through solo-LTR formation affected mostly older elements. On the other hand, Ivana, very much like Tork, shows a low number of solo-LTRs and a high number of young intact elements, suggesting active retrotransposition and low recombination rates. In contrast, Angela, Sire, TAR, Bianca and Ikeros families have more solo-LTRs than intact elements, although the numbers are all low, suggesting limited retrotransposition activity for each of these families. Finally, Ale, Galadriel, and Reina show the same number of intact LTR-RTs and solo-LTRs in the tetraploids as in the parental diploids, and similar LTR-RT age, suggesting little or no activity during the recent evolution of \u003cem\u003eGossypium\u003c/em\u003e.\u003c/p\u003e\n\u003ch3\u003eCRM and Tekay insertion and elimination have shaped the centromeres of tetraploid cottons\u003c/h3\u003e\n\u003cp\u003eThe CRM and Tekay lineages have previously been shown to have amplified in the tetraploids [\u003cspan citationid=\"CR26\" class=\"CitationRef\"\u003e26\u003c/span\u003e] and have been suggested to play a role in reorganization of their centromeres [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e]. Our results confirm that CRM and Tekay have been actively transposing after the polyploidization event, as deduced by their recent insertion time distribution (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e). Moreover, our results also show that CRM copy number increased greatly in the A subgenomes of both tetraploids (up to 8-fold in \u003cem\u003eG. hirsutum\u003c/em\u003e and up to 6.5-fold in \u003cem\u003eG. barbadense\u003c/em\u003e) and to a lesser extent in the D subgenomes (up to 2.6-fold in \u003cem\u003eG. hirsutum\u003c/em\u003e and up to 2-fold in \u003cem\u003eG. barbadense\u003c/em\u003e), with this increase being slightly higher in \u003cem\u003eG. hirsutum\u003c/em\u003e than in \u003cem\u003eG. barbadense\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e top left). A different pattern was observed for Tekay elements, which exhibited increased copy numbers in the D subgenomes of both tetraploids (also slightly more in \u003cem\u003eG. hirsutum\u003c/em\u003e than in \u003cem\u003eG. barbadense\u003c/em\u003e), but not in the A subgenomes, where its copy number remained high and stable in \u003cem\u003eG. hirsutum\u003c/em\u003e and decreased about 30% in \u003cem\u003eG. barbadense\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e top right). These results thus show that, whereas both CRM and Tekay actively transposed accompanying polyploidization, their amplification rate was dissimilar and differed between species and subgenomes.\u003c/p\u003e \u003cp\u003eAnalysis of the CRM and Tekay solo-LTRs shows that the elimination rate also varies between the two polyploids and between subgenomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e bottom panel). The number of CRM solo-LTRs is much higher in the D subgenomes than in the A subgenomes in both polyploids, but it is particularly high in \u003cem\u003eG. barbadense\u003c/em\u003e. This could be the consequence of a higher elimination rate of recent CRM insertions from the D subgenomes, particularly in \u003cem\u003eG. barbadense\u003c/em\u003e, which would also explain the lower number of intact CRMs in the D subgenomes of this species (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e top left) and the slightly older insertion time of the remaining \u003cem\u003eG. barbadense\u003c/em\u003e CRMs (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e middle panel). Interestingly, the analysis of the Tekay lineage suggests an opposite dynamic. Tekay solo-LTRs are much more abundant in the A subgenomes of both tetraploids than in the D subgenomes, a trend that is even more extreme for \u003cem\u003eG. hirsutum\u003c/em\u003e (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e bottom right). Interestingly, the high number of Tekay solo-LTRs, doubling that of complete elements, is also seen in the diploid parental, \u003cem\u003eG. herbaceum\u003c/em\u003e, which may suggest a link between the nature of the A subgenomes and the efficiency of Tekay recombination. This high elimination rate of Tekay elements in the A subgenomes could explain why the recent activity suggested by the very recent age of the intact elements (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e middle panel) has only resulted in an increase of Tekay insertion in the D subgenomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig3\" class=\"InternalRef\"\u003e3\u003c/span\u003e top panel).\u003c/p\u003e \u003cp\u003eThe distribution of the intact CRM elements and solo-LTR CRM shows a high concentration of CRM intact copies in the regions characterized as containing the centromeres [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e], whereas CRM solo-LTRs are more homogeneously distributed along the chromosomes (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003ea, Supplementary Fig.\u0026nbsp;2), suggesting that, although CRM elements may target the centromere regions [\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e], their differential elimination from other chromosome regions reinforces their concentration in the centromeric regions. Similarly, intact Tekay elements in the D subgenomes are concentrated in the pericentromeric regions whereas Tekay solo-LTRs are not, suggesting active elimination of Tekay elements from chromosomal regions maintains their prevalence in pericentromeric regions. Interestingly, Tekay elements are not found in CRM-rich regions of the centromere, perhaps reflecting an insertional preference difference between these two types of elements. The distribution of Tekay intact and solo-LTRs in the A subgenomes seems less skewed than in the D subgenome, although they also show an opposite distribution with respect to genes.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eComparison of the CRM and Tekay distributions between the diploids and the tetraploids also reveals notable differences. First, the high increase in CRM copy number in the A subgenomes of the two polyploids, and particularly in \u003cem\u003eG. hirsutum\u003c/em\u003e (Supplementary Fig.\u0026nbsp;3), is concentrated almost exclusively in the centromere, whereas in the diploids, the distribution of the CRM elements is less skewed (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eb, and Supplementary Fig.\u0026nbsp;3). Second, the solo-LTR distribution in the diploids parallels that of the intact copies, suggesting that there is no preferential elimination of CRM elements from chromosome arms (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eb). For Tekay elements, the increase in intact elements in the D subgenomes of polyploids (up to 5-fold in G. hirsutum and up to 3.5-fold in G. barbadense) without a concomitant increase in number of solo-LTRs appears to concentrate the intact elements in the pericentromeric regions (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eb).\u003c/p\u003e\n\u003ch3\u003eA small number of LTR-RT lineages accounts for intraspecies polymorphisms\u003c/h3\u003e\n\u003cp\u003eAlthough the number of insertions of most LTR-RT lineages appears relatively constant within each tetraploid species, the young age of the insertions for many LTR-RT lineages suggests recent retrotransposition. Therefore, to understand the scale and scope of TE activity within species, we constructed reference-based pangenomes for both species that we subsequently used to characterize transposon insertion polymorphisms (TIPs) within \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e. We used genome assemblies of seven different varieties to build the \u003cem\u003eG.hirsutum\u003c/em\u003e pangenome and genome assemblies of nine varieties and nine different varieties to build the \u003cem\u003eG. barbadense\u003c/em\u003e pangenome. The different genomes sequences were aligned to the references for each species (TM-1 and 3\u0026ndash;79) to build each pangenome. The alignments covered 94.1\u0026ndash;98.7% the \u003cem\u003eG. hirsutum\u003c/em\u003e TM1 reference and 94.8\u0026ndash;97.0% of the \u003cem\u003eG. barbadense\u003c/em\u003e 3\u0026ndash;79 reference. The number of structural variants (SV) detected in the \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e pangenomes was 63,081 and 55,122, respectively, most of which (57,727 and 48,990, respectively) were insertion/deletion (indel) polymorphisms. Approximately 20% of these indels (12,316 \u003cem\u003ein G. hirsutum\u003c/em\u003e and 9,386 in \u003cem\u003eG. barbadense\u003c/em\u003e) contained LTR-RT sequences (at least 30% coverage of an LTR-RT). To detect recent LTR-RT insertions, we defined LTR-RTs TIPs as indel SVs with at least 80% identity over 80% of the length of an intact LTR-RT representative or a solo-LTR, and a minimum size of 80% of an intact LTR-RT or a solo-LTR. This analysis led to the detection of 4436 LTR-RT TIPs in \u003cem\u003eG. hirsutum\u003c/em\u003e and 3520 in \u003cem\u003eG. barbadense\u003c/em\u003e, and 2406 solo-LTR TIPs in \u003cem\u003eG. hirsutum\u003c/em\u003e and 1313 in \u003cem\u003eG. barbadense\u003c/em\u003e. These SVs were homogeneously distributed along chromosomes, as were the LTR-RTs and solo-LTR TIPs (Supplementary Fig.\u0026nbsp;4a).\u003c/p\u003e \u003cp\u003eIn addition to this overall quantitative view of TE-mediated TIPS, LTR-RT lineage classification of both intact and solo-LTR TIPs shows that only few of the many LTR-RT lineages are responsible for the polymorphism generated within the two cotton species (Figs.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea and \u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eb, Supplementary Tables\u0026nbsp;1,2). Interestingly, some lineages that experienced copy number increase in tetraploid cotton relative to their diploid progenitors, e.g., Ogre (Supplementary Fig.\u0026nbsp;1), also show little variability within the polyploids (less than 1%, Supplementary Tables\u0026nbsp;1,2 and Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea), suggesting that they experienced a burst of transposition following polyploidization but became silent before the split of the two tetraploids from their common ancestor. In contrast, some lineages appear to have retained high transpositional activity after speciation of the two polyploids. These include the \u003cem\u003eGypsy\u003c/em\u003e CRM and Tekay lineages and the \u003cem\u003eCopia\u003c/em\u003e Ivana and Tork lineages (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003ea). Interestingly, CRM is mainly polymorphic as solo-LTRs, whereas Ivana and Tork are mainly polymorphic as complete LTR-RTs and Tekay is notably polymorphic as both. These patterns could suggest recent elimination activity for CRM, recent retrotransposition for Ivana and Tork, and insertion followed by elimination of new copies for Tekay. Whereas CRM and Tekay elements concentrate in centromeric and pericentromeric regions (see above), \u003cem\u003eCopia\u003c/em\u003e elements Ivana and Tork have a more widespread chromosomal distribution (not shown) and are more likely to have been inserted near genes and regulatory regions. We explored the potential effect of these polymorphic Tork and Ivana insertions by evaluating their position relative to genes in the pangenomes of both species (Fig.\u0026nbsp;\u003cspan refid=\"Fig5\" class=\"InternalRef\"\u003e5\u003c/span\u003eb). Our analysis shows that Ivana and Tork insertions are significantly enriched inside and near (2 kbp) genes in both species and may therefore have impacted gene coding capacity and regulation, potentially triggering phenotypic consequences.\u003c/p\u003e \u003cp\u003e \u003c/p\u003e\n\u003ch3\u003eLTR-RT TIPs impact on cotton domestication and breeding\u003c/h3\u003e\n\u003cp\u003eTo address the possibility that genomic variability triggered by recent TIPs had an impact on cotton domestication and improvement, we analyzed their presence in wild and domesticated \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e populations, including landraces and elite varieties. To this end we mapped short-read resequencing data of 283 varieties of \u003cem\u003eG. hirsutum\u003c/em\u003e and 223 varieties of \u003cem\u003eG. barbadense\u003c/em\u003e to their respective pangenome graphs and genotyped the presence/absence of the 6841 \u003cem\u003eG. hirsutum\u003c/em\u003e and 5749 \u003cem\u003eG. barbadense\u003c/em\u003e TIPs across these populations. The varieties analyzed represent the wild-to-domesticated continuum of both species including, 151 and 51 elite varieties, 66 and 60 landraces, and 73 and 112 wild accessions, of \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e, respectively (Supplementary Table\u0026nbsp;3) [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e].\u003c/p\u003e \u003cp\u003ePrincipal Component Analysis (PCA) shows that TIPs provide enough signal to distinguish between groups of wild, landrace, and domesticated varieties and suggests that there has been extensive differential transposition activity and/or retention/elimination of LTR-RT copies and solo-LTRs between the three groups of accessions (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ea,b). This is particularly clear for \u003cem\u003eG. hirsutum\u003c/em\u003e, in line with previous analyses based on SNPs which showed that the cultivated \u003cem\u003eG. hirsutum\u003c/em\u003e accessions clustered tightly due to their narrow genetic diversity [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e]. This drastic reduction of diversity accompanying domestication was not seen in \u003cem\u003eG. barbadense\u003c/em\u003e, probably because modern \u003cem\u003eG. barbadense\u003c/em\u003e cultivars have a more complex and obscure origin that in addition involves many intentional introgressions [\u003cspan citationid=\"CR24\" class=\"CitationRef\"\u003e24\u003c/span\u003e].\u003c/p\u003e \u003cp\u003e \u003c/p\u003e \u003cp\u003eAn analysis of the TIP frequency in G. \u003cem\u003ehirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e populations shows that there are few LTR-RT insertions or deletions present in only one or two groups (see Supplementary Table\u0026nbsp;4) whereas most TIPs (5860 in \u003cem\u003eG. hirsutum\u003c/em\u003e and 4450 in \u003cem\u003eG. barbadense\u003c/em\u003e) are present in all three groups (i.e., wild, landrace, and domesticated). Interestingly, however, the TIPs shared by the three groups are present in very different frequencies among accession in those groups (Supplementary Figs.\u0026nbsp;5 and 6). The presence of TIPs at different population frequencies suggests the possibility of selection, positive or purifying, acting on polymorphic LTR-RTs or on closely linked \u003cem\u003eloci\u003c/em\u003e. To find signs indicative of positive selection of TIPs in cultivated varieties, we used the population branch statistics (PBS) method. This approach measures the level of genetic differentiation for a specific population along a phylogenetic branch, providing insights into population-specific evolutionary changes [\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e]. Figure\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ec and \u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ed show the PBS distribution of LTR-RT insertions and deletions relative to the reference genome in cultivated versus landrace and wild accessions and in wild versus cultivated and landrace accessions in \u003cem\u003eG. hirsutum.\u003c/em\u003e\u003c/p\u003e \u003cp\u003eWe identified the TIPs with the highest PBS (above 95th percentile in the comparisons cultivated against landraces and wild, and wild against cultivated and landraces). Eighty percent of the high-PBS TIPs of the wild versus cultivated and landraces (linked to domestication) and 63% of high-PBS TIPs in cultivated versus wild and landraces (linked to cotton breeding) are related to \u003cem\u003eGypsy\u003c/em\u003e LTR-RT (mainly Tekay) elements (Supplementary Table\u0026nbsp;5). Tekay elements are concentrated in pericentromeric region of chromosomes, particularly in the D subgenomes, and our results show that they have contributed to the recent evolution of these regions, as discussed above. Many of the Tekay TIPs with high PBS are far from genes (Supplementary Table\u0026nbsp;4), have a distribution similar to the Tekay intact elements (Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003e and Supplementary Fig.\u0026nbsp;4b), and are in regions characterized by a low density in coding genes. The results presented here suggest that Tekay insertions in these regions may have been selected in the recent evolution of \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e. Notably, although most high-PBS Tekay TIPs are far from genes, others are located close to genes (less than 2 Kb away, Supplementary Table\u0026nbsp;5 and Suplementary Fig.\u0026nbsp;\u003cspan refid=\"Fig4\" class=\"InternalRef\"\u003e4\u003c/span\u003eb) and may have impacted their coding capacity or regulation. This is also the case of High-PBS TIPs related to \u003cem\u003eCopia\u003c/em\u003e LTR-RTs, mainly the Ivana and Tork superfamilies, which are frequently found tightly associated with genes.\u003c/p\u003e \u003cp\u003eThe high PBS associated with these TIPs may suggest selection due to favorable or detrimental effects on the nearby genes. We analyzed the expression of genes tightly linked with high-PBS TIPs using previously published data for the \u003cem\u003eG. hirsutum\u003c/em\u003e TM1 accession [\u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e] and found that they show different patterns of expression throughout cotton development and under different stress conditions (Supplementary Fig.\u0026nbsp;7). An analysis of the functions of these genes shows a significant enrichment for \"auxin response\" genes (p-value\u0026thinsp;=\u0026thinsp;0,00039, GO term GO:0009733) in \u003cem\u003eG. hirsutum\u003c/em\u003e. We found 5 genes belonging to this functional category (four coding for a SAUR-like auxin-responsive proteins and one for an Indoleacetic acid-induced protein 16) within the 95 percentile of highest PBS in \u003cem\u003eG. hirsutum\u003c/em\u003e, indicative of high population frequency differentiation. This represents an 11-fold enrichment in comparison to the total annotated predicted genes (Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ee). Two of these TIPs are shown in Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eg and \u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003eh. For the first of these, we observed a strong population frequency differentiation of the landraces and cultivated accessions with respect to the wild (potential relationship with domestication), whereas for the second, we found a strong differentiation between the cultivated and the other two groups, suggesting that this TIP could have been targeted during crop improvement. A preliminary analysis of the TIPs with high PBS value in \u003cem\u003eG. barbadense\u003c/em\u003e TIPs also showed an enrichment for the function \"response to auxin\u0026rdquo; (p\u0026thinsp;=\u0026thinsp;0.022, Fig.\u0026nbsp;\u003cspan refid=\"Fig6\" class=\"InternalRef\"\u003e6\u003c/span\u003ee, Supplementary Table\u0026nbsp;4), suggesting that variants in auxin-related genes were selected in the domestication and breeding process of the two species of cotton. Whether this result reflects parallel selection or historical interspecific introgression is at present an open question.\u003c/p\u003e"},{"header":"Discussion","content":"\u003cp\u003ePlant TEs, and in particular LTR-RTs, have been traditionally assumed to evolve through bursts of amplification followed by periods of low activity and elimination [\u003cspan citationid=\"CR11\" class=\"CitationRef\"\u003e11\u003c/span\u003e]. Recent data, however, suggest that the diversity of LTR-RTs that coexist within plant genomes often results from heterogeneous evolutionary dynamics [\u003cspan citationid=\"CR12\" class=\"CitationRef\"\u003e12\u003c/span\u003e]. The results presented here show that LTR-RT elements are in general more abundant and recent in the two tetraploid cotton species compared with models of their diploid progenitors, suggesting a general burst of amplification accompanying the polyploidization event, as it is often but not always the case in plants [\u003cspan citationid=\"CR9\" class=\"CitationRef\"\u003e9\u003c/span\u003e]. However, these results also show that different LTR-RT lineages experienced rather different evolutionary dynamics during cotton polyploidization, speciation, and subsequent crop improvement processes. As an example, whereas Ogre elements amplified in the tetraploids, they probably became silent before speciation, as they show little polymorphism between the two tetraploid species. Other LTR-RT lineages, such as Tork and Ivana, have probably retained activity for longer, since more than half of their insertions are polymorphic within the species. This is also the case for CRM and Tekay, which have recently been proposed to play important roles in the evolution of \u003cem\u003eGossypium\u003c/em\u003e centromeres [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e, \u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eOur results show that CRM elements amplified after polyploidization. This CRM element number increase has been higher in the A than the D subgenomes, probably because they are more efficiently removed from the former subgenome, as evidenced by a higher number of solo-LTRs. As CRMs concentrate in cotton centromeric regions, this leads to a differential increase of CRMs in the centromeric regions of the A subgenomes in the two tetraploid species as compared to the model diploid progenitor genome (represented by modern \u003cem\u003eG. herbaceum\u003c/em\u003e). Moreover, the number of CRM solo-LTRs in the D subgenomes is higher in \u003cem\u003eG. barbadense\u003c/em\u003e than in \u003cem\u003eG. hirsutum\u003c/em\u003e, which correlates with a lower number of intact CRMs in this species, suggesting that CRM elimination from the centromeric regions of the D subgenomes is higher in this species. In contrast to CRM LTR-RTs, Tekay elements have greatly increased in the D subgenomes while their number has remained constant in the A subgenomes. These results are in line with recent evidence regarding the role of these two lineages in centromere evolution in cotton after polyploidization, with the proportion of CRM versus Tekay elements being higher in the centromeric regions of the A subgenomes, and lower in those of the D subgenomes [\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e]. Our results add to these earlier observations in suggesting that this shift may reflect differences in elimination of the insertions between the two subgenomes rather than differential insertion. Moreover, the analysis of CRM and Tekay intact and solo-LTR distributions along chromosomes suggests that, although these elements may be centrophilic[\u003cspan citationid=\"CR32\" class=\"CitationRef\"\u003e32\u003c/span\u003e] and target the centromere for integration, their elimination via solo-LTR formation from chromosome arms reinforces their concentration in centromeric and pericentromeric regions. Taken together, our results are compatible with the important role of CRM and Tekay LTR-RTs in the expansion and reorganization of cotton centromeres after polyploidization that has been recently suggested[\u003cspan citationid=\"CR27\" class=\"CitationRef\"\u003e27\u003c/span\u003e, \u003cspan citationid=\"CR28\" class=\"CitationRef\"\u003e28\u003c/span\u003e] and points to the role of their selective elimination from chromosomal regions and subgenomes in specifying their dynamics and their genomic distribution.\u003c/p\u003e \u003cp\u003eWith respect to continuing LTR-RT activity \u003cem\u003ewithin\u003c/em\u003e species, our results show that CRM and especially Tekay intact LTR-RT and solo-LTR insertions are highly polymorphic within the two polyploid species, and are present at different population frequencies in wild, landrace and cultivar groups of accessions. Interestingly, many of the Tekay TIPs are among those with the highest Population Branch Statistic (PBS) values, suggesting potential selection during domestication and breeding. The potential selection of Tekay TIPs sitting in pericentromeric regions suggests that variants of these heterochromatic regions could have selectable phenotypic impacts. Centromeres are fast evolving structures[\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e] and it has previously been shown that centromere variants have been selected during the domestication of crops, such as maize and wheat [\u003cspan citationid=\"CR33\" class=\"CitationRef\"\u003e33\u003c/span\u003e, \u003cspan citationid=\"CR34\" class=\"CitationRef\"\u003e34\u003c/span\u003e]. The selection of centromere variants could be both due to differences in the centromere function or due to the selection of particular alleles of genes residing in this low recombining region. Clarifying this point will require further analyses, but the results presented here point to an important role of genetic variability in heterochromatic regions for crop domestication and breeding.\u003c/p\u003e \u003cp\u003eIn addition to being important players in centromere evolution [\u003cspan citationid=\"CR29\" class=\"CitationRef\"\u003e29\u003c/span\u003e], LTR-RTs have also played a major role as drivers of genetic diversity in plant gene coding capacity and regulation [\u003cspan citationid=\"CR2\" class=\"CitationRef\"\u003e2\u003c/span\u003e, \u003cspan citationid=\"CR35\" class=\"CitationRef\"\u003e35\u003c/span\u003e]. Our data show that, in addition to \u003cem\u003eGypsy\u003c/em\u003e CRM and Tekay, two \u003cem\u003eCopia\u003c/em\u003e LTR-RT lineages (Tork and Ivana) have amplified and remained active in the two polyploid species of cotton since their divergence from a common ancestor. Notably, we identified a large number of very recent retrotranspositions, as evidenced by their identical LTRs, that are significantly enriched within exons and in the proximal regions upstream and downstream of genes. Many of these insertions have varying frequencies in the wild, landrace, and cultivar groups, exhibiting high PBS values and suggesting that selection during cotton domestication and crop improvement may have shaped their distribution. Gene expression data reveals a wide range of developmental and stress-related expression patterns. Considering the many different traits targeted by domestication and crop improvement, including plant architecture, fruiting habit, flower and seed development and fiber-length and quality [\u003cspan citationid=\"CR36\" class=\"CitationRef\"\u003e36\u003c/span\u003e], more work is needed to determine whether the genomic changes introduced these genic and gene-adjacent LTR-RT insertions have contributed to these processes. It is interesting to note that both \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e exhibited an overrepresentation of auxin-responsive genes, and in particular SAUR-like genes, among those impacted by LTR-RT insertions. Auxins are well-established as key controllers of plant development and stress responses [\u003cspan citationid=\"CR37\" class=\"CitationRef\"\u003e37\u003c/span\u003e], which also play a role in the regulation cotton fiber development [\u003cspan citationid=\"CR38\" class=\"CitationRef\"\u003e38\u003c/span\u003e, \u003cspan citationid=\"CR39\" class=\"CitationRef\"\u003e39\u003c/span\u003e] and broader processes such as plant growth [\u003cspan citationid=\"CR40\" class=\"CitationRef\"\u003e40\u003c/span\u003e, \u003cspan citationid=\"CR41\" class=\"CitationRef\"\u003e41\u003c/span\u003e]. Moreover, SAUR genes, which are key players in plant adaptation and growth [\u003cspan citationid=\"CR42\" class=\"CitationRef\"\u003e42\u003c/span\u003e] have been shown to have been targets of domestication and breeding of crops such as rice, citrus, and \u003cem\u003eBrassica oleracea\u003c/em\u003e [\u003cspan additionalcitationids=\"CR44\" citationid=\"CR43\" class=\"CitationRef\"\u003e43\u003c/span\u003e\u0026ndash;\u003cspan citationid=\"CR45\" class=\"CitationRef\"\u003e45\u003c/span\u003e]. The potential independent selection of transposon insertion polymorphisms closely linked to different auxin-related genes in both \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e, highlights the importance of this hormone signaling in cotton domestication and breeding and suggests convergent selection in both species (notwithstanding the possibility of an interspecific introgressive origin).\u003c/p\u003e"},{"header":"Conclusion","content":"\u003cp\u003eThis study highlights the importance of considering lineage level TE classifications when studying LTR-RT dynamics, as well as examining not only retrotransposition but also recombination and the formation of solo-LTRs as the joint proximal determinants of TE presence in genomes. Our study shows LTR-RT dynamics are highly differential among LTR-RT families, and that only some were activated after polyploidisation and disparately between the two cotton species and the two subgenomes in each. CRM and Tekay elements appear to have played key roles in centromere reorganization, and Tork and Ivana are suggested to actively generate variability within and close to genes. These few recently active lineages have generated LTR-RT polymorphisms that may have been selected during cotton domestication and crop improvement, have modified chromosome heterochromatic regions and have impacted genes expressed in different tissues and environmental situations, including genes related to auxin signaling, which are key players in plant and fiber development regulation, in both \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e.\u003c/p\u003e"},{"header":"Methods","content":"\u003cp\u003eAll scripts and pipelines used for this study are available within the following directory: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://github.com/Lcamdom/panTEvo\u003c/span\u003e\u003cspan address=\"https://github.com/Lcamdom/panTEvo\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/p\u003e \u003cdiv id=\"Sec11\" class=\"Section2\"\u003e \u003ch2\u003eAnnotation of intact LTR-RTs and solo LTRs\u003c/h2\u003e \u003cp\u003eAll genomes selected and included in this study (supplementary table 6) were assembled into pseudo-chromosomes. These were annotated for intact LTR-RT elements using EDTA-raw scripts on the LTR-mode [\u003cspan citationid=\"CR46\" class=\"CitationRef\"\u003e46\u003c/span\u003e]. The intact elements found were further classified using TEsorter [\u003cspan citationid=\"CR47\" class=\"CitationRef\"\u003e47\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eAll intact elements within the same species and lineage were clustered using CD-hit [\u003cspan citationid=\"CR48\" class=\"CitationRef\"\u003e48\u003c/span\u003e] for 80% similarity and 80% coverage. The representatives for each cluster were used to build species-specific, home-made in silico libraries of intact LTR-RT elements. For each element in these libraries, the left LTR was extracted using the structural annotation generated by LTR_Retriever [\u003cspan citationid=\"CR49\" class=\"CitationRef\"\u003e49\u003c/span\u003e], and species-specific home-made libraries of LTRs were generated. Solo-LTRs were then identified in each genome using these libraries. Blast [\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e] was used to find LTRs across the genome, and solo-LTRs were then defined by filtering out any hit with any LTR-RT sequence 1kb either upstream or downstream.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec12\" class=\"Section2\"\u003e \u003ch2\u003ePangenome construction and TIP characterisation\u003c/h2\u003e \u003cp\u003eMinimap2 [\u003cspan citationid=\"CR51\" class=\"CitationRef\"\u003e51\u003c/span\u003e] was used to map every \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e genome to the references (TM1 and 3\u0026ndash;79, respectively). Structural variants (SVs) were called using svim-asm [\u003cspan citationid=\"CR52\" class=\"CitationRef\"\u003e52\u003c/span\u003e] and the vcfs were merged using bcftools [\u003cspan citationid=\"CR53\" class=\"CitationRef\"\u003e53\u003c/span\u003e] and collapsed with truvari [\u003cspan citationid=\"CR54\" class=\"CitationRef\"\u003e54\u003c/span\u003e]. The pangenome graphs were generated using the vg toolkit [\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eTIPs were characterised using blast [\u003cspan citationid=\"CR50\" class=\"CitationRef\"\u003e50\u003c/span\u003e] to find homology between the pangenome SVs and the intact LTR-RT and solo-LTR libraries. Resulting hits were excluded when they represented less than 80% of the length of the TE consensus, had lower than 80% homology with the reference TE, and covered less than 80% of the SV.\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec13\" class=\"Section2\"\u003e \u003ch2\u003ePangenome genotyping and population frequency analysis\u003c/h2\u003e \u003cp\u003eRe-sequencing data from 295 and 222 accessions across the wild-to-domesticated continuum (Supplementary Table\u0026nbsp;3) were mapped to the \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e pangenomes, respectively. Mapping of reads to the pangenome graph and structural variant calling were accomplished using the vg toolkit [\u003cspan citationid=\"CR55\" class=\"CitationRef\"\u003e55\u003c/span\u003e].\u003c/p\u003e \u003cp\u003eThe TIP population frequencies were calculated using R. TIPs present in the three population groups were combined in a single matrix to obtain PBS values. TIPs with MAF\u0026thinsp;\u0026gt;\u0026thinsp;0.01 and call rate\u0026thinsp;\u0026gt;\u0026thinsp;95% were used to calculate Fst values using SNPready [\u003cspan citationid=\"CR56\" class=\"CitationRef\"\u003e56\u003c/span\u003e]. PBS values were calculated following the formulas described in [\u003cspan citationid=\"CR30\" class=\"CitationRef\"\u003e30\u003c/span\u003e]. GO enrichment of genes with proximal TIPs with high PBS was performed with GOATOOLS [\u003cspan citationid=\"CR57\" class=\"CitationRef\"\u003e57\u003c/span\u003e].\u003c/p\u003e \u003c/div\u003e \u003cdiv id=\"Sec14\" class=\"Section2\"\u003e \u003ch2\u003eExpression analysis of genes hear high PBS TIPs\u003c/h2\u003e \u003cp\u003e \u003cem\u003eG. hirsutum\u003c/em\u003e RNA-seq data [\u003cspan citationid=\"CR31\" class=\"CitationRef\"\u003e31\u003c/span\u003e] was processed to estimate transcript abundance using Trinity Transcript Quantification scripts [\u003cspan citationid=\"CR58\" class=\"CitationRef\"\u003e58\u003c/span\u003e] and Salmon [\u003cspan citationid=\"CR59\" class=\"CitationRef\"\u003e59\u003c/span\u003e] was used as the abundance estimation method. The raw count matrices were processed into Rlog values using DESeq2 [\u003cspan citationid=\"CR60\" class=\"CitationRef\"\u003e60\u003c/span\u003e]. Heatmaps are normalised per row (scale\u0026thinsp;=\u0026thinsp;rows), and were generated using the Pheatmap function in R.\u003c/p\u003e \u003c/div\u003e"},{"header":"Declarations","content":"\u003ch2\u003eAuthor Contribution\u003c/h2\u003e\u003cp\u003eLCD and JMC designed the project with the help of all other authors. LCD did most of the analyses with the help of RC for TIP genotyping and TIP population analyses. JFW and CEG contributed cotton sequence information and resequencing data for cotton varieties. LCD and JMC wrote the paper with contributions ofrom all other authors.\u003c/p\u003e\u003ch2\u003eData Availability\u003c/h2\u003e\u003cp\u003eThe source code supporting this study is available in the GitHub repository https://github.com/Lcamdom/panTEvo.\u003c/p\u003e"},{"header":"References","content":"\u003col\u003e\u003cli\u003e\u003cspan\u003eWendel JF, Jackson SA, Meyers BC, Wing RA. Evolution of plant genome architecture. Genome Biol [Internet]. 2016;17:37. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://genomebiology.biomedcentral.com/articles/10.1186/s13059-016-0908-1\u003c/span\u003e\u003cspan address=\"http://genomebiology.biomedcentral.com/articles/10.1186/s13059-016-0908-1\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLisch D. How important are transposons for plant evolution? Nat Rev Genet [Internet]. 2013;14:49\u0026ndash;61. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.ncbi.nlm.nih.gov/pubmed/23247435\u003c/span\u003e\u003cspan address=\"http://www.ncbi.nlm.nih.gov/pubmed/23247435\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAndersson L, Purugganan M. Molecular genetic variation of animals and plants under domestication. Proc Natl Acad Sci U S A. 2022;119:e2122150119.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLiu P, Cuerda-Gil D, Shahid S, Slotkin R. The Epigenetic Control of the Transposable Element Life Cycle in Plant Genomes and Beyond. Annu Rev Genet. 2022;56:63\u0026ndash;87.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLuo S, Zhang H, Duan Y, Yao X, Clark AG, Lu J. The evolutionary arms race between transposable elements and piRNAs in Drosophila melanogaster. BMC Evol Biol [Internet]. 2020;20:14. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1186/s12862-020-1580-3\u003c/span\u003e\u003cspan address=\"10.1186/s12862-020-1580-3\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKosuge M, Ito J, Hamada M. Landscape of evolutionary arms races between transposable elements and KRAB-ZFP family. Sci Rep [Internet]. 2024;14:23358. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1038/s41598-024-73752-7\u003c/span\u003e\u003cspan address=\"10.1038/s41598-024-73752-7\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLawlor MA, Ellison CE. Evolutionary dynamics between transposable elements and their host genomes: mechanisms of suppression and escape. Curr Opin Genet Dev [Internet]. 2023;82:102092. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.sciencedirect.com/science/article/pii/S0959437X23000722\u003c/span\u003e\u003cspan address=\"https://www.sciencedirect.com/science/article/pii/S0959437X23000722\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang H, Lang Z, Zhu JK. Dynamics and function of DNA methylation in plants. Nat Rev Mol Cell Biol. 2018. pp. 489\u0026ndash;506.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVicient CM, Casacuberta JM. Impact of transposable elements on polyploid plant genomes. Ann Bot [Internet]. 2017;120:195\u0026ndash;207. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://dx.doi.org/10.1093/aob/mcx078\u003c/span\u003e\u003cspan address=\"10.1093/aob/mcx078\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNieto Feliner G, Casacuberta J, Wendel JF. Genomics of Evolutionary Novelty in Hybrids and Polyploids. Front Genet. 2020;11:792.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBaidouri M, El, Panaud O. Comparative genomic paleontology across plant kingdom reveals the dynamics of TE-driven genome evolution. Genome Biol Evol [Internet]. 2013;5:954\u0026ndash;65. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3673626\u0026amp;tool=pmcentrez\u0026amp;rendertype=abstract\u003c/span\u003e\u003cspan address=\"http://www.pubmedcentral.nih.gov/articlerender.fcgi?artid=3673626\u0026amp;tool=pmcentrez\u0026amp;rendertype=abstract\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePulido M, Casacuberta JM. Transposable element evolution in plant genome ecosystems. Curr Opin Plant Biol. 2023;75:102418.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStritt C, Thieme M, Roulin AC. Rare transposable elements challenge the prevailing view of transposition dynamics in plants. Am J Bot. 2021;108:1310\u0026ndash;4.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStitzer MC, Anderson SN, Springer NMV, Ross-Ibarra J. The Genomic Ecosystem of Transposable Elements in Maize. PLoS Genet. 2021;17:e1009768.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePereira V. Insertion bias and purifying selection of retrotransposons in the Arabidopsis thaliana genome. Genome Biol. 2004;5:R79.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePiegu B, Guyot R, Picault N, Roulin A, Saniyal A, Kim H, et al. Doubling genome size without polyploidization: Dynamics of retrotransposition-driven genomic expansions in Oryza australiensis, a wild relative of rice. Genome Res. 2006;16:1262\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBennetzen JL, Kellogg EA. Do Plants Have a One-Way Ticket to Genomic Obesity? Plant Cell [Internet]. 1997;9:1509\u0026ndash;14. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.plantcell.org/content/9/9/1509.short\u003c/span\u003e\u003cspan address=\"http://www.plantcell.org/content/9/9/1509.short\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMunasinghe M, Read A, Stitzer M, Song B, Menard C, Ma K, et al. Combined analysis of transposable elements and structural variation in maize genomes reveals genome contraction outpaces expansion. PLoS Genet. 2023;19:e1011086.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDevos KM, Brown JKM, Bennetzen JL. Genome size reduction through illegitimate recombination counteracts genome expansion in Arabidopsis. Genome Res. 2002;12:1075\u0026ndash;9.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTian Z, Rizzon C, Du J, Zhu L, Bennetzen JL, Jackson SA, et al. Do genetic recombination and gene density shape the pattern of DNA elimination in rice long terminal repeat retrotransposons? Genome Res. 2009;19:2221\u0026ndash;30.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShirasu K, Schulman A, Lahaye T, Schulze-Lefert P. A contiguous 66-kb barley DNA sequence provides evidence for reversible genome expansion. Genome Res. 2000;10:908\u0026ndash;15.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eVitte C, Panaud O. Formation of solo-LTRs through unequal homologous recombination counterbalances amplifications of LTR retrotransposons in rice Oryza sativa L. Mol Biol Evol. 2003;20(4):528\u0026thinsp;\u0026ndash;\u0026thinsp;40. Mol Biol Evol. 2003;20:528\u0026ndash;40.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eViot CR, Wendel JF. Evolution of the Cotton Genus, Gossypium, and Its Domestication in the Americas. CRC Crit Rev Plant Sci [Internet]. 2023;42:1\u0026ndash;33. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1080/07352689.2022.2156061\u003c/span\u003e\u003cspan address=\"10.1080/07352689.2022.2156061\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDaojun Y, Grover CE, Hu G, Pan M, Miller ER, Conover JL, et al. Parallel and Intertwining Threads of Domestication in Allopolyploid Cotton. Adv Sci. 2021;8:2003634.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang M, Li J, Wang P, Liu F, Liu Z, Zhao G, et al. Comparative Genome Analyses Highlight Transposon-Mediated Genome Expansion and the Evolutionary Architecture of 3D Genomic Folding in Cotton. Mol Biol Evol. 2021;38:3621\u0026ndash;36.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHe X, Qi Z, Liu Z, Chang X, Zhang X, Li J et al. Pangenome analysis reveals transposon-driven genome evolution in cotton. BMC Biol [Internet]. 2024;22:92. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1186/s12915-024-01893-2\u003c/span\u003e\u003cspan address=\"10.1186/s12915-024-01893-2\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChang X, He X, Li J, Liu Z, Pi R, Luo X, et al. High-quality Gossypium hirsutum and Gossypium barbadense genome assemblies reveal the landscape and evolution of centromeres. Plant Commun. 2024;5:100722.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi J, Liu Z, You C, Qi Z, You J, Grover CE et al. Convergence and divergence of diploid and tetraploid cotton genomes. Nat Genet. 2024;10.1038/s41588-024-01964\u0026ndash;8.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eNaish M, Henderson IR. The structure, function, and evolution of plant centromeres. Genome Res. 2024;34:161\u0026ndash;78.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eYi X, Liang Y, Huerta-Sanchez E, Jin1 X, Cuo Z, Ping X et al. Sequencing of Fifty Human Exomes Reveals Adaptation to High Altitude. Science (1979). 2010;329:75\u0026ndash;78.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHu Y, Chen J, Fang L, Zhang Z, Ma W, Niu Y et al. Gossypium barbadense and Gossypium hirsutum genomes provide insights into the origin and evolution of allotetraploid cotton. Nature Genetics. 2019 51:4 [Internet]. 2019 [cited 2025 Feb 21];51:739\u0026ndash;48. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.nature.com/articles/s41588-019-0371-5\u003c/span\u003e\u003cspan address=\"https://www.nature.com/articles/s41588-019-0371-5\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eBousios A, Kakutani T, Henderson IR. Centrophilic Retrotransposons of Plant Genomes. Annu Rev Plant Biol [Internet]. 2025 [cited 2025 Mar 6]; Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/39952673/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/39952673/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eSchneider KL, Xie Z, Wolfgruber TK, Presting GG. Inbreeding drives maize centromere evolution. Proc Natl Acad Sci U S A [Internet]. 2016 [cited 2025 Feb 14];113:E987\u0026ndash;96. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/26858403/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/26858403/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eWang Z, Wang W, Xie X, Wang Y, Yang Z, Peng H et al. Dispersed emergence and protracted domestication of polyploid wheat uncovered by mosaic ancestral haploblock inference. Nat Commun [Internet]. 2022 [cited 2025 Feb 14];13:3891. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pmc.ncbi.nlm.nih.gov/articles/PMC9259585/\u003c/span\u003e\u003cspan address=\"https://pmc.ncbi.nlm.nih.gov/articles/PMC9259585/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMeyer RS, Purugganan MD. Evolution of crop species: Genetics of domestication and diversification. Nat Rev Genet [Internet]. 2013;14:840\u0026ndash;52. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttp://www.ncbi.nlm.nih.gov/pubmed/24240513\u003c/span\u003e\u003cspan address=\"http://www.ncbi.nlm.nih.gov/pubmed/24240513\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGrover CE, Yoo MJ, Lin M, Murphy MD, Harker DB, Byers RL et al. Genetic Analysis of the Transition from Wild to Domesticated Cotton (Gossypium hirsutum L.). G3: Genes|Genomes|Genetics [Internet]. 2019 [cited 2025 Feb 21];10:731. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pmc.ncbi.nlm.nih.gov/articles/PMC7003101/\u003c/span\u003e\u003cspan address=\"https://pmc.ncbi.nlm.nih.gov/articles/PMC7003101/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eTeale WD, Paponov IA, Palme K. Auxin in action: signalling, transport and the control of plant growth and development. Nature Reviews Molecular Cell Biology 2006 7:11 [Internet]. 2006 [cited 2025 Feb 21];7:847\u0026ndash;59. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.nature.com/articles/nrm2020\u003c/span\u003e\u003cspan address=\"https://www.nature.com/articles/nrm2020\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eJareczek JJ, Grover CE, Wendel JF. Cotton fiber as a model for understanding shifts in cell development under domestication. Front Plant Sci [Internet]. 2023;14. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.frontiersin.org/journals/plant-science/articles/\u003c/span\u003e\u003cspan address=\"https://www.frontiersin.org/journals/plant-science/articles/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003e10.3389/fpls.2023.1146802\u003c/span\u003e\u003cspan address=\"10.3389/fpls.2023.1146802\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eChu Q, Fu X, Zhao J, Li Y, Liu L, Zhang L et al. Simultaneous improvement of fiber yield and quality in upland cotton (Gossypium hirsutum L.) by integration of auxin transport and synthesis. Molecular Breeding [Internet]. 2024 [cited 2025 Feb 21];44:1\u0026ndash;17. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://link.springer.com/article/10.1007/s11032-024-01500-w\u003c/span\u003e\u003cspan address=\"https://link.springer.com/article/10.1007/s11032-024-01500-w\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eMa J, Pei W, Ma Q, Geng Y, Liu G, Liu J et al. QTL analysis and candidate gene identification for plant height in cotton based on an interspecific backcross inbred line population of Gossypium hirsutum \u0026times; Gossypium barbadense. Theoretical and Applied Genetics [Internet]. 2019 [cited 2025 Feb 21];132:2663\u0026ndash;76. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://link.springer.com/article/\u003c/span\u003e\u003cspan address=\"https://link.springer.com/article/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003e10.1007/s00122-019-03380-7\u003c/span\u003e\u003cspan address=\"10.1007/s00122-019-03380-7\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eShi H, Qanmber G, Yang Z, Guo Y, Ma S, Shu S, et al. An AP2/ERF transcription factor GhERF109 negatively regulates plant growth and development in cotton. Plant Sci. 2025;352:112365.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eStortenbeker N, Bemer M. The SAUR gene family: the plant\u0026rsquo;s toolbox for adaptation of growth and development. J Exp Bot [Internet]. 2019 [cited 2025 Feb 24];70:17\u0026ndash;27. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/30239806/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/30239806/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGonzalez-Ibeas D, Ibanez V, Perez-Roman E, Borred\u0026aacute; C, Terol J, Talon M. Shaping the biology of citrus: II. Genomic determinants of domestication. Plant Genome [Internet]. 2021 [cited 2025 Feb 24];14. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/34464512/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/34464512/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhao Z, Chen T, Yue J, Pu N, Liu J, Luo L et al. Small Auxin Up RNA 56 (SAUR56) regulates heading date in rice. Mol Breed [Internet]. 2023 [cited 2025 Feb 24];43. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/37521314/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/37521314/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGuo N, Wang S, Wang T, Duan M, Zong M, Miao L, et al. A graph-based pan-genome of Brassica oleracea provides new insights into its domestication and morphotype diversification. Plant Commun. 2024;5:100791.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOu S, Su W, Liao Y, Chougule K, Agda JRA, Hellinga AJ, et al. Benchmarking transposable element annotation methods for creation of a streamlined, comprehensive pipeline. Genome Biol. 2019;20:275.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eZhang R-G, Li G-Y, Wang X-L, Dainat J, Wang Z-X, Ou S et al. TEsorter: An accurate and fast method to classify LTR-retrotransposons in plant genomes. Hortic Res [Internet]. 2022;9:uhac017. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1093/hr/uhac017\u003c/span\u003e\u003cspan address=\"10.1093/hr/uhac017\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eFu L, Niu B, Zhu Z, Wu S, Li W. CD-HIT: Accelerated for clustering the next-generation sequencing data. Bioinformatics. 2012;28:3150\u0026ndash;2.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eOu S, Jiang N, LTR_retriever:. A Highly Accurate and Sensitive Program for Identification of Long Terminal Repeat Retrotransposons. Plant Physiol [Internet]. 2018 [cited 2025 Jan 7];176:1410\u0026ndash;22. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/29233850/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/29233850/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eAltschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403\u0026ndash;10.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLi H. Minimap2: Pairwise alignment for nucleotide sequences. Bioinformatics. 2018;34:3094\u0026ndash;100.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHeller D, Vingron M. SVIM-asm: structural variant detection from haploid and diploid genome assemblies. Bioinformatics [Internet]. 2021;36:5519\u0026ndash;21. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1093/bioinformatics/btaa1034\u003c/span\u003e\u003cspan address=\"10.1093/bioinformatics/btaa1034\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eDanecek P, Bonfield JK, Liddle J, Marshall J, Ohan V, Pollard MO et al. Twelve years of SAMtools and BCFtools. Gigascience [Internet]. 2021;10:giab008. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1093/gigascience/giab008\u003c/span\u003e\u003cspan address=\"10.1093/gigascience/giab008\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eEnglish AC, Menon VK, Gibbs RA, Metcalf GA, Sedlazeck FJ. Truvari: refined structural variant comparison preserves allelic diversity. Genome Biol [Internet]. 2022;23:271. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1186/s13059-022-02840-6\u003c/span\u003e\u003cspan address=\"10.1186/s13059-022-02840-6\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eHickey G, Heller D, Monlong J, Sibbesen JA, Sir\u0026eacute;n J, Eizenga J et al. Genotyping structural variants in pangenome graphs using the vg toolkit. Genome Biol [Internet]. 2020;21:35. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://doi.org/10.1186/s13059-020-1941-7\u003c/span\u003e\u003cspan address=\"10.1186/s13059-020-1941-7\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGranato ISC, Galli G, de Oliveira Couto EG, e Souza MB, Mendon\u0026ccedil;a LF, Fritsche-Neto R. snpReady: a tool to assist breeders in genomic analysis. Molecular Breeding [Internet]. 2018 [cited 2025 Jan 7];38:1\u0026ndash;7. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://link.springer.com/article/\u003c/span\u003e\u003cspan address=\"https://link.springer.com/article/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003e10.1007/s11032-018-0844-8\u003c/span\u003e\u003cspan address=\"10.1007/s11032-018-0844-8\" targettype=\"DOI\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eKlopfenstein DV, Zhang L, Pedersen BS, Ram\u0026iacute;rez F, Vesztrocy AW, Naldi A et al. GOATOOLS: A Python library for Gene Ontology analyses. Scientific Reports 2018 8:1 [Internet]. 2018 [cited 2025 Jan 7];8:1\u0026ndash;17. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://www.nature.com/articles/s41598-018-28948-z\u003c/span\u003e\u003cspan address=\"https://www.nature.com/articles/s41598-018-28948-z\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eGrabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, et al. Full-length transcriptome assembly from RNA-Seq data without a reference genome. Nat Biotechnol. 2011;29:644\u0026ndash;52.\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003ePatro R, Duggal G, Love MI, Irizarry RA, Kingsford C. Salmon provides fast and bias-aware quantification of transcript expression. Nat Methods [Internet]. 2017 [cited 2025 Feb 27];14:417\u0026ndash;9. Available from: \u003cspan class=\"ExternalRef\"\u003e\u003cspan class=\"RefSource\"\u003ehttps://pubmed.ncbi.nlm.nih.gov/28263959/\u003c/span\u003e\u003cspan address=\"https://pubmed.ncbi.nlm.nih.gov/28263959/\" targettype=\"URL\" class=\"RefTarget\"\u003e\u003c/span\u003e\u003c/span\u003e\u003c/span\u003e\u003c/li\u003e \u003cli\u003e\u003cspan\u003eLove MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550.\u003c/span\u003e\u003c/li\u003e\u003c/ol\u003e"}],"fulltextSource":"","fullText":"","funders":[],"hasAdminPriorityOnWorkflow":false,"hasManuscriptDocX":true,"hasOptedInToPreprint":true,"hasPassedJournalQc":"","hasAnyPriority":true,"hideJournal":false,"highlight":"","institution":"","isAcceptedByJournal":true,"isAuthorSuppliedPdf":false,"isDeskRejected":"","isHiddenFromSearch":false,"isInQc":false,"isInWorkflow":false,"isPdf":false,"isPdfUpToDate":true,"isWithdrawnOrRetracted":false,"journal":{"display":true,"email":"[email protected]","identity":"genome-biology","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"gbio","sideBox":"Learn more about [Genome Biology](https://genomebiology.biomedcentral.com/)","snPcode":"13059","submissionUrl":"https://submission.springernature.com/new-submission/13059/3","title":"Genome Biology","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"stoa","reportingPortfolio":"BMC/SO AJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true},"keywords":"retrotransposon, solo-LTR, allopolyploid plants, pangenome, TIP (Transposon Insertion Polymorphism), evolutionary genomics","lastPublishedDoi":"10.21203/rs.3.rs-6172192/v1","lastPublishedDoiUrl":"https://doi.org/10.21203/rs.3.rs-6172192/v1","license":{"name":"CC BY 4.0","url":"https://creativecommons.org/licenses/by/4.0/"},"manuscriptAbstract":"\u003ch2\u003eBackground\u003c/h2\u003e \u003cp\u003eTransposable elements (TEs) are major components of plant genomes and major drivers of plant genome evolution. The cotton genus (\u003cem\u003eGossypium\u003c/em\u003e) is an excellent evolutionary model for polyploidization, speciation, domestication and crop improvement. Here, we implement genome and pangenome analyses to study in detail the dynamics of LTR-retrotransposons (LTR-RT) during the cotton evolution.\u003c/p\u003e\u003ch2\u003eResults\u003c/h2\u003e \u003cp\u003eWe show that some LTR-RT lineages amplified in tetraploid cotton compared to their diploid progenitors, whereas others stayed stable or amplified but were removed through solo-LTR formation. Using species-level pangenomes we show that only a few lineages (CRM, Tekay, Ivana and Tork) remained active after polyploidization and are still transposing. Tekay and CRM elements have re-shaped the centromeric and pericentromeric regions of tetraploid cottons in a subgenome specific manner, through new insertions but also selective eliminations through solo-LTR formation. On the other hand, Ivana and Tork have actively inserted within or close to genes. Finally, population-level analyses using the two pangenomes and data from 283 and 223 varieties of \u003cem\u003eG. hirsutum\u003c/em\u003e and \u003cem\u003eG. barbandense\u003c/em\u003e reveal changes in Transposon Insertion Polymorphism (TIP) frequencies accompanying domestication and improvement of both species, suggesting the possibility of selection on linked regions.\u003c/p\u003e\u003ch2\u003eConclusions\u003c/h2\u003e \u003cp\u003eOur findings reveal that LTR-RT lineages followed differential dynamics during cotton evolution, displaying differences among species and the two coresident genomes of allopolyploid cotton. A handful of the LTR-RT lineages that expanded after polyploidisation helped shape the genomes of both \u003cem\u003eG. hirstutum\u003c/em\u003e and \u003cem\u003eG. barbadense\u003c/em\u003e, impacting their centromere and pericentromeric regions as well as protein- coding genes.\u003c/p\u003e","manuscriptTitle":"Differential LTR-retrotransposon dynamics across polyploidization, speciation, domestication and improvement of cotton (Gossypium)","msid":"","msnumber":"","nonDraftVersions":[{"code":1,"date":"2025-03-21 05:34:30","doi":"10.21203/rs.3.rs-6172192/v1","editorialEvents":[{"type":"communityComments","content":0},{"type":"decision","content":"Revision requested","date":"2025-05-21T13:17:58+00:00","index":"","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-05-05T15:11:08+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-05-02T17:55:30+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-04-26T11:12:31+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"333169007904495181083168643754359692955","date":"2025-04-18T08:10:17+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"91871230115627795878906315024705479555","date":"2025-04-16T02:37:19+00:00","index":"hide","fulltext":""},{"type":"editorInvitedReview","content":"","date":"2025-04-11T15:19:45+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"44274741331874411523976744621244431525","date":"2025-04-01T17:45:36+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"216307206451589781546238925592001614386","date":"2025-04-01T16:51:48+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"13407148536033180419071140657502403629","date":"2025-03-31T19:42:24+00:00","index":"hide","fulltext":""},{"type":"reviewerAgreed","content":"211707649843810873909734619868891015125","date":"2025-03-31T15:35:28+00:00","index":"hide","fulltext":""},{"type":"reviewersInvited","content":"","date":"2025-03-31T15:03:48+00:00","index":"","fulltext":""},{"type":"editorAssigned","content":"","date":"2025-03-17T11:43:49+00:00","index":"","fulltext":""},{"type":"checksComplete","content":"","date":"2025-03-10T10:36:44+00:00","index":"","fulltext":""},{"type":"submitted","content":"Genome Biology","date":"2025-03-06T15:57:10+00:00","index":"","fulltext":""}],"status":"published","journal":{"display":true,"email":"[email protected]","identity":"genome-biology","isNatureJournal":false,"hasQc":true,"allowDirectSubmit":false,"externalIdentity":"gbio","sideBox":"Learn more about [Genome Biology](https://genomebiology.biomedcentral.com/)","snPcode":"13059","submissionUrl":"https://submission.springernature.com/new-submission/13059/3","title":"Genome Biology","twitterHandle":"","acdcEnabled":true,"dfaEnabled":true,"editorialSystem":"stoa","reportingPortfolio":"BMC/SO AJ","inReviewEnabled":true,"inReviewRevisionsEnabled":true}}],"origin":"","ownerIdentity":"a982337e-2302-44a0-a98b-0e9f17ee2a52","owner":[],"postedDate":"March 21st, 2025","published":true,"recentEditorialEvents":[],"rejectedJournal":[],"revision":"","amendment":"","status":"published-in-journal","subjectAreas":[],"tags":[],"updatedAt":"2025-11-03T16:05:02+00:00","versionOfRecord":{"articleIdentity":"rs-6172192","link":"https://doi.org/10.1186/s13059-025-03837-7","journal":{"identity":"genome-biology","isVorOnly":false,"title":"Genome Biology"},"publishedOn":"2025-10-27 15:58:08","publishedOnDateReadable":"October 27th, 2025"},"versionCreatedAt":"2025-03-21 05:34:30","video":"","vorDoi":"10.1186/s13059-025-03837-7","vorDoiUrl":"https://doi.org/10.1186/s13059-025-03837-7","workflowStages":[]},"version":"v1","identity":"rs-6172192","journalConfig":"researchsquare"},"__N_SSP":true},"page":"/article/[identity]/[[...version]]","query":{"redirect":"/article/rs-6172192","identity":"rs-6172192","version":["v1"]},"buildId":"XKTyCvWXoU3ODBz1xrDgd","isFallback":false,"isExperimentalCompile":false,"dynamicIds":[84888],"gssp":true,"scriptLoader":[]}