HelitronScanner uncovers a large overlooked cache of Helitron transposons in many plant genomes

Helitron Terminal Features Represented by LCVs.

Based on an improved LCV algorithm (22), HelitronScanner aims to extract more definitive Helitron features than the few previously identified: the TC dinucleotide at the 5′ end, the hairpin structure and CTRR (R = A or G) sequence at the 3′ end, and the A and T residues flanking the 5′ and 3′ ends, respectively. More than 5,000 double-ended Helitron sequences in maize, Arabidopsis, rice, sorghum, Caenorhabditis elegans, and Medicago truncatula genomes were taken to create two sets of 100-bp slices from both Helitron ends, which were then clustered separately by using cd-hit (24) at 90% similarity to remove redundant sequences and, thus, avoid LCV bias. A training set was created from 2,846 seed sequences of 5′-end clusters and 2,048 of 3′-end clusters to extract LCVs (Methods), and 1,613 double-ended Helitrons were assembled from these paired seed sequences to evaluate LCVs (see Confidence Level by LCV Scores below). LCVs from Helitron termini convey local conserved information in the sense that no global multiple sequence alignment is required that, in turn, overcomes the drawback of overall Helitron variation and gives rise to a more generalized set of features for Helitron identification in a broader array of organisms. LCVs combined together as identification features reveal more terminal conservation than previous methods that focus on overall sequence alignment.

Sequence logos of 303 LCVs from the Helitron 5′-terminal 50 bp (Fig. 1A) and 575 LCVs from the 3′-terminal 50 bp (Fig. 1B) display patterns in sequence conservation, compared with sequence logos from the 5,676 raw sequences of the same regions (Fig. 1 C and D). Helitron 5′ ends appear to be AT-rich, whereas 3′ ends contain a higher content of Cs and Gs between 11 and 30 bp from the endpoint. It is interesting that the dinucleotide AA at the 21 and 20 positions in Fig. 1B coincides with the dinucleotide AA often found in the middle of the CG base-paired hairpin loops (14) and that the hairpin structure at the 3′ end is reflected in LCVs without any prior knowledge input. Because of significant sequence diversity among Helitrons, there is little conserved sequence information other than 5′-TC and 3′-CTAG (Fig. 1 C and D).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 1.

Helitron terminal conservation shown in sequence logos. Sequence logos of: 303 LCVs generated from 50 bp at Helitron 5′ termini (A); 575 LCVs generated from 50 bp at Helitron 3′ termini (B); raw sequences of Helitron 5′ termini (C); raw sequences of Helitron 3′ termini (D). LCV location variation (A and B) has been normalized according to their frequency.

The LCV algorithm extracts overrepresented patterns based on an iteration of exhaustive enumeration of oligonucleotides. The patterns may have nonconserved locations and do not have to be vertically aligned to a fixed location, which brings great flexibility for motif finding at the cost of speed. After trials of LCV extraction, we determined from the location of generated LCVs that the Helitron’s conserved regions are 50 bp at 5′ and 3′ ends, so LCVs were extracted from them. The 10 most frequent LCVs from 5,676 Helitrons are shown in Table 1. LCV notation follows the syntax of regular expression in computer science, i.e., a dot denotes any nucleotide, brackets denote alternative nucleotides, and numbers within braces denote occurrences. For instance, AA.G.ACG.{9}CT[AG]{2} represents a motif comprised of 2 As, a nucleotide of any kind (i.e., N), 1 G, 1 N, ACG, 9 Ns, and the CTRR (R for A or G) ending. An LCV with nonconserved regions of fixed length, represented by a single number within braces, provides more discriminative information than an LCV with nonconserved regions of variable length. Among the 5,676 Helitrons in the training set, the frequencies of 3′ LCVs are higher than those of 5′ LCVs, indicating that Helitron 3′ ends are more conserved than 5′ ends, in agreement with our previous study (14).

Confidence Level by LCV Scores.

The sum of the LCV scores of the predicted Helitron’s 5′ and 3′ ends is taken as the primary criterion of the prediction’s confidence level. Higher scores indicate more specificity of the combination of LCV features and, thus, provide a higher confidence level. However, too stringent a threshold may cause more diverse Helitrons to be missed. After analyzing the distribution of LCV scores among the 1,613 Helitrons in the training set, a threshold of 5 for the LCV score at each end was chosen initially as a compromise between sensitivity and specificity (SI Appendix, Fig. S1). Statistically, 95% of the Helitrons in the training set had scores ≥10.

Helitron Identification in Angiosperms.

Using an LCV score threshold of 5 for each end, we ran HelitronScanner against a wide range of plant genome sequences from Phytozome version 9.0 (23) and identified 107,367 Helitrons. The LCV scores assigned by HelitronScanner to the identified Helitrons are an indicator of prediction confidence. SI Appendix, Fig. S2 shows the variation of predicted Helitron number from all plant genomes under different LCV score thresholds. The number of Helitrons decreases dramatically with more stringent LCV score thresholds. There are 107,367 Helitrons using an LCV score ≥10 as the cutoff criterion, 33,530 (or 31.2% of original) at ≥20, 12,439 (11.6%) at ≥30, 7,812 (7.3%) at ≥40, and 5,164 (4.8%) at ≥50. Higher LCV scores provide a higher prediction confidence, whereas more stringent thresholds lead to the loss of a fraction of true Helitrons (see In Silico Verification of Helitrons below).

SI Appendix, Table S1 shows the number of Helitrons in each organism, their genome abundance, and size distribution. Because HelitronScanner identifies Helitron 5′ and 3′ termini separately and the 3′ termini are more conserved, the Helitron number is mainly determined by the 3′ termini, whereas the 5′ termini are more diverse and abundant. For each species, we then estimated a false positive rate (FPR) at the chosen LCV threshold by running HelitronScanner on a randomized version of its genome (Methods). The FPR varied from species to species, and several of those lacking close relatives in the training set often had an estimated FPR >50%. A more stringent 3′-end threshold score of 10 was chosen for those species to balance detection power against FPR (SI Appendix, Table S2). Plant species, such as algae and bryophytes, still giving an estimated FPR >45% under the new threshold, were dropped from the list. Table 2 gives the number of Helitrons predicted by HelitronScanner and their percentage in angiosperm genomes, after correcting for estimated false positives.

A thorough search for Helitrons in the maize genome by HelitronFinder (12) and HelSearch (11) identified 2,791 and 1,930 Helitrons, respectively, which comprised approximately 2% of the genome. In contrast, HelitronScanner identified 31,233 Helitrons or 6.6% of the maize genome, including 92% of those identified by HelitronFinder. To verify new Helitrons uniquely identified by HelitronScanner, we carried out tests of insertion site polymorphisms by PCR assays and in silico comparisons. The results indicate that HelitronScanner is efficacious in identifying authentic Helitrons (see details in Helitron verification). HelSearch (11) identified 281 complete Helitrons in A. thaliana, 230 in Medicago, 651 in rice, and 608 in sorghum, compared with HelitronScanner, which identified 609, 1,142, 4,634, and 2,082 Helitrons in the respective genomes. Three rice accessions (japonica, indica, and glaberrima) were used to test Helitron insertion polymorphisms: 248 of the Helitrons predicted by HelitronScanner were verified, in contrast to 179 of those predicted by HelSearch.

HelitronScanner also identified Helitrons in many plant species where Helitrons had not been previously reported (organisms marked with an asterisk in Table 2), including a wide range of monocots and eudicots. The Helitron abundance revealed by our study is consistent with previous suggestions that the percentage of Helitrons in genomes was most likely underestimated (12). We found no sign of correlation between Helitron abundance and genome sizes.

Helitron Verification by PCR Assays in Multiple Maize Lines.

Maize has the most variable genome structure yet described (25, 26) and one of the highest contents of transposons (85%) among fully sequenced genomes, components of the genome that have shaped its architecture over time (27). Different maize inbreds are highly polymorphic in their transposon content and distribution and provide valuable germplasm for the validation of computationally predicted Helitrons (14). Often, a Helitron is present in some lines while absent in others. Here, we also adopted the plus/minus polymorphism criterion to validate the authenticity of predicted Helitrons.

We randomly picked 15 high-score (LCV > 20) and 4 medium-score (LCV = 11–20) Helitrons predicted exclusively by HelitronScanner (i.e., not by HelitronFinder or HelSearch) that were flanked by single-copy regions in the B73 reference genome (28). PCR primers for flanking sequences were designed and plus/minus variation was tested in different inbred lines (SI Appendix, Table S3 from hel_pcr_001 to hel_pcr_019). Of the 19 Helitrons, 13 exhibited vacant sites (i.e., only flanking sequences amplified) in maize lines other than B73, thus verifying their authenticity (Column “Validated”; TRUE entries in SI Appendix, Table S3). A PCR band image of five validated Helitrons is also shown in SI Appendix, Fig. S3. Helitrons not validated here by the PCR assay are not necessarily false positives because they may be absent in inbred lines not included in our small panel.

A supporting criterion of a Helitron’s authenticity is element copy number. As seen in SI Appendix, Table S3, either a high LCV score or a high copy number (>4 copies in the host genome) is a strong predictor of Helitron authenticity.

In Silico Verification of Helitrons.

Applying the same concept as in the PCR verification, we also compared Helitron insertion sites and their flanking sequences in silico among multiple sequenced accessions where true Helitrons might exhibit plus/minus polymorphism. By BLASTing Helitron flanking sequence joints against other sequenced accessions, we verified Helitrons in silico based on the presence of vacant sites. SI Appendix, Fig. S4 shows examples of Helitron polymorphisms detected in maize, rice, and Arabidopsis. As can be seen, the flanking regions are highly conserved, so vacant sites lacking Helitron insertions are easily identified. We tested maize Helitrons predicted in the fully assembled genome of the inbred B73 (28) against contigs of the Mo17 inbred (http://bo.csam.montclair.edu/du/software/helitronscanner). SI Appendix, Fig. S4A shows that a 1,572-bp Helitron identified on chromosome 1 of B73 was absent in Mo17. In rice, the fully assembled genomes of the two subspecies japonica (29) and indica (30) of Oryza sativa and contig data of Oryza glaberrima from Arizona Genomics Institute were used to test Helitron polymorphism. SI Appendix, Fig. S4B shows that a 2,009-bp Helitron identified on chromosome 1 of japonica was absent in glaberrima. In Arabidopsis, the ecotypes Columbia, C24, and Bur-0 from the A. thaliana 1001 Genomes Project were investigated. As seen in SI Appendix, Fig. S4C, a 547-bp Helitron was present in Col and C24, but not in Bur-0. The identified Helitron and its flanking region are highly conserved in Col and C24. Compared with other programs, HelitronScanner identified many inactive Helitrons that are too divergent to have been detected previously, which, in turn, caused a decrease in the Helitron verification rate. For example, HelSearch detected 651 full-length Helitrons and at least 6,947 elements with conserved 3′ ends in the rice subspecies japonica. The detection power of HelitronFinder, however, is highly confined to the maize genome.

To validate Helitrons efficiently, artificial sequences of 50 bp were made up by joining 25-bp flanking sequences on both sides of the Helitron and BLASTed against genome sequences for evidence of plus or minus polymorphisms. The mega-BLAST task of nucleotide blastn was chosen to search for highly similar sequences, and a minimum 80% coverage of the query sequence was required to support the presence of the joint, i.e., vacant, flanking sequences. Within one organism, the presence of vacant sequences in accessions other than the one with the predicted Helitron validated the Helitron’s authenticity. We compared numbers of validated Helitrons identified by HelitronScanner with those identified by HelitronFinder or HelSearch, two widely used computer programs, in maize, rice, and Arabidopsis (Table 3). A large number of maize Helitrons were missed by the previous methods, showing the efficacy of HelitronScanner in identifying Helitrons in the highly polymorphic maize genome. In rice, HelitronScanner also predicted 69 more validated Helitrons than HelSearch. However, in Arabidopsis, HelitronScanner identified three fewer validated Helitrons than HelSearch simply because these Helitrons had 3′ LCV scores below 5, the HelitronScanner threshold chosen to avoid a high false positive rate.

We analyzed the LCV scores of 1,616 Helitrons validated in maize, rice and Arabidopsis to assess retrospectively our selection of a threshold (SI Appendix, Table S4). Because Helitrons scoring lower than 10 were not included, based on the distribution of scores in the training set and PCR assays, the minimum score is 10. Of the validated Helitrons, only 15 (0.9%) have an LCV score of 10, 152 (9.4%) have a score of 11–20, and the vast majority (89.7%) have scores >20 (Fig. 2). The maximum, average, and median LCV scores are 75, 48, and 56, respectively. In particular, 813 (50.3%) of the validated Helitrons have scores higher than 50. The score distribution indicates that validated Helitrons mainly have high scores and that low-scoring Helitrons are rare. Copy number is a complementary indicator when Helitron scores are low, because most Helitrons have at least two copies. Of the validated Helitrons, only 30 are singleton and their scores are all >10. Therefore, a pragmatic guideline would be to accept Helitrons with LCV scores >20 and multiple copies, reject singleton Helitrons with LCV scores ≤10, and analyze intermediate cases further.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 2.

Score distribution of 1,616 validated Helitrons.

Evolutionary Distance Revealed by LCVs.

LCVs are overrepresented patterns attributed to Helitrons. We analyzed the more conserved 3′ termini by showing how LCVs are shared among species in terms of evolutionary distance (Fig. 3). We looked for the presence of the extracted 575 LCVs from the 3′ ends in the 1,616 in silico-validated Helitrons (1,352 from maize, 248 from rice, and 16 from Arabidopsis). A 1,616 × 575 matrix was generated based on the matching condition “1 for true and 0 for false” of every Helitron against every LCV, which was then decomposed by principal component analysis. Each Helitron was projected onto the top-two principal components (PC) 1 and PC2, which accounted for 15.03% and 3.07% of total sample variance, respectively (Fig. 3). Although Helitrons in the monocots (maize and rice) overlap, those in the eudicot Arabidopsis are distinctive. This interrelationship reveals that LCVs cannot only serve collectively as the deterministic feature of Helitrons, but can also convey evolutionary relationships.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 3.

Clustering of validated Helitrons by LCV principal components. Five hundred seventy-five LCVs from the 3′ ends of Helitrons in the training set were matched against 1,616 in silico-validated Helitrons, including 1,352 maize Helitrons (red circles) 248 rice Helitrons (blue circles), and 16 Helitrons from A. thaliana (green triangles). All Helitrons were projected to the top two significant principal components, PC1 and PC2, which account for 15.03% and 3.07% of total sample variance, respectively.

More interestingly, the pattern GC.CG.{9}CTRR is the most shared sequence feature among all studied species. In particular, 959 (70.9%) of 1,352 validated (i.e., polymorphic and, therefore, recently transposed) maize Helitrons and 156 (62.9%) of 248 validated rice Helitrons contain this pattern, which may reflect a crucial characteristic for Helitron transposition.