De novo, high-quality assembly and annotation of the halophyte grass Aeluropus littoralis draft genome and identification of A20/AN1 zinc finger protein family

Abstract

Background

Aeluropus littoralis is considered a valuable natural forage plant for ruminant livestock and is highly tolerant to extreme abiotic stresses, especially salinity, drought, and heat. It is a monocotyledonous halophyte, has salt glands, performs C4-type photosynthesis and has a close genetic relationship with cereal crops. Moreover, previous studies have shown its huge potential as a reservoir of genes and promoters to understand and improve abiotic stress tolerance in crops.

Results

The sequencing and hybrid assembly of the A. littoralis genome (2n = 2X = 20) using short and long reads from the BGISeq-500 and PacBio high-fidelity (HiFi) sequencing platforms, respectively. Using the k-mer analysis method, the haploid genome size of A. littoralis was estimated to be 360 Mb (with a heterozygosity rate of 1.88%). The hybrid assembled genome included 4,078 contigs with a GC content of 44% and covered 348 Mb. The longest contig and the N50 values were 5.1 Mb and 133.77 kb, respectively. The Benchmarking Universal Single Copy Ortholog (BUSCO) value was 91.1%, indicating good integrity of the assembled genome. The discovered repetitive elements accounted for 90.6 Mb, representing 26.03% of the total genome, and included a significant component of transposable elements (11.48%, ~40 Mb). Using a homology-based approach, 35,147 genes were predicted from the genome assembly. We next focused our analysis on the zinc-finger A20/AN1 gene family, a member of which (AlSAP) was previously shown to confer increased tolerance to osmotic and salt stresses when it was over-expressed in tobacco, wheat, and rice. Here, we identified the complete set of members of this family in the Aeluropus littoralis genome, thereby laying the foundation for their future functional analysis in cereal crops. In addition, the expression patterns of four novel genes from this family were analyzed by qPCR.

Conclusion

This resource and our findings will contribute to improve the current understanding of salinity tolerance in halophytes while providing useful genes and allelic variation to improve salinity and drought tolerance in cereals through genetic engineering and gene editing, respectively.

Introduction

Salinization has emerged as a critical issue for agricultural propagation as it affects 20% of irrigated lands globally particularly in arid and semi-arid regions [1,2,3]. Plants are constantly exposed to various biotic and abiotic stresses, which greatly impact their growth and productivity [4, 5]. To cope with these threats, plants trigger genetic networks of molecular pathways that promote stress tolerance [6,7,8,9,10]. These adaptive responses involve dynamic adjustments in morphology, metabolism, and physiology throughout the plant's life cycle [11, 12]. Over the past two decades, significant progress has been made in elucidating the transcriptional changes induced by abiotic stresses, as well as identifying the signaling proteins and transcription factors that regulate stress-responsive gene expression [13]. In the Mediterranean and Middle East, water scarcity and soil salinity pose significant challenges to food security, especially cereal production. To address these challenges, modern biotechnological approaches, including molecular markers, genetic engineering, genomics, molecular physiology, and whole-genome sequencing, are being utilized to enhance the tolerance of cereals like wheat, barley, and rice to salt, drought, and heat stress [14]. A promising strategy within this framework involves harnessing the genetic traits of naturally resilient plants, such as halophytes and xerophytes, which have evolved mechanisms to thrive in saline and arid environments [15, 16]. Halophyte plant species are characterized by their remarkable ability to grow and multiply in highly saline conditions (≥200 mM NaCl) and are thus preferred over glycophytes [15, 17, 18]. These plants exhibit various morphological adaptations to counter abiotic stress, particularly salinity, including reduced leaf surface area to limit transpiration, leaf rolling, salt secretion via specialized glands or bladders, and deep, extensive root systems that access lower-salinity water sources [19]. Consequently, halophytes are considered sustainable, low-cost candidates for improving soil physical and chemical properties [20, 21].

Aeluropus littoralis, a perennial halophytic monocot grass from the Poaceae family, has recently attracted attention owing to its utility as valuable natural forage and its remarkable tolerance to environmental stresses [22,23,24]. Taxonomically, it belongs to the subtribe Aeluropodinae, within the tribe Cynodonteae and the subfamily Chloridoideae, a member of the PACMAD clade, which includes the majority of C4 photosynthetic plants [25]. Phylogenetic analyses based on chloroplast genomes have shown that A. littoralis, A. lagopoides, and A. sinensis form a single subtribe, Aeluropodinae, which is a sister to the Euleusininae subtribe. Notably, A. littoralis is closely related to Eleusine coracana and Eleusine indica (subtribe Euleusininae) [25]. Ecologically, this species thrives in salt marshes and exhibits high tolerance to both drought- and heat stress [22, 26]. Furthermore, A. littoralis has also been identified as a valuable reservoir of stress-responsive genes and promoters, with potential applications in enhancing abiotic stress resistance in commercial crops [26,27,28,29]. High-throughput screening of suppression subtractive hybridization (SSH) libraries prepared from salt-stressed roots and leaves of A. littoralis identified 492 differentially expressed sequence tags (ESTs) [22]. Among them, the AlSAP gene, encoding a zinc finger A20/AN1 stress-associated protein (SAP), has been functionally characterized. Its overexpression in transgenic tobacco, rice, and durum wheat has led to enhanced tolerance to salt and osmotic stress through regulation of stress-related gene networks [26, 27, 30,31,32,33].

Recent advancements in genome sequencing technologies and assembly algorithms have significantly improved the construction of high-quality, chromosome-scale plant genomes. The integration of long-read sequencing platforms (e.g., PacBio, Oxford Nanopore) with proximity-based methods like Hi-C or classical genetic mapping has enabled the precise scaffolding of sequence contigs, thereby improving genome assembly enhancing accuracy and completeness [34]. In contrast, short-read assemblers are not efficient for resolving complex and repetitive regions such as centromeres and telomeres [35]. Hybrid assembly strategies that integrate long read (Pacific Biosciences (PacBio)/Oxford Nanopore) and short read (Illumina/BGI) sequencing data have demonstrated superior effectiveness in assembling more contiguous plant genomes [34,35,36]. To date, over 1,139 plant genomes have been published using such techniques [37]. The genomes of several orphan grasses with notable tolerance to drought and salinity have recently been reported, including Eleusine coracana [38], the halophyte Oryza coarctata [39], Eleusine indica [40], Eragrostis curvula [41], Puccinellia tenuiflora [42], and Achnatherum splendens [43]. For A. littoralis, only a single draft de novo genome assembly (Iranian ecotype) has been published, based on short reads from the Illumina HiSeq 2500 platform [24]. This assembly was highly fragmented, comprising 182,747 contigs with a contig N50 of just 3.6 kb and covering 300 Mb of the estimated 354 Mb genome, as determined by flow cytometry [44]. Moreover, it included only 15,916 predicted genes models-considerably fewer than in related monocots- while repetitive elements comprised just 21.6% of the genome, significantly lower than the average repeat content (45.49%) observed across plant genomes [45].

In the present study, we analyzed the chromosome number and estimated the genome size of A. littoralis (Tunisian ecotype) using k-mer analysis method. Notably, we generated a high-quality draft genome assembly of this diploid halophytic C4 grass using a hybrid strategy combining of long reads (PacBio) and short reads (BGISeq-500). The genome assembly integrity was evaluated through BUSCO analysis and by mapping 431 of the 492 previously identified salt-related ESTs [22]. The resulting genome assembly enabled the identification and characterization of SAP gene family members encoding A20/AN1 zinc finger proteins. This genomic resource provides a valuable foundation for understanding the molecular basis of abiotic stress tolerance in A. littoralis, and offers promising tools for improving stress resistance in commercial crops via gene transfer and CRISPR/Cas-based genome editing.

Results

Morphological and anatomical responses of A. littoralis to salt stress

As shown in Fig. 1, increasing NaCl concentrations resulted in significant reductions in the shoot length, leaf size, and internode length by up to 60% at 550 mM NaCl relative to plants grown under normal conditions (Fig. S1). Additionally, root system architecture was markedly altered, with a noticeable decrease in both lateral and crown roots (CRs) numbers. Interestingly, root length increased by approximately 70 % at moderate salinity levels (150–300 mM NaCl) relative to the control (Fig. 1A). Conversely, at concentrations exceeding 300 mM NaCl, roots growth declined by 30–35%, aligning more closely with glycophytic rather than halophytic responses (Fig. S1), A noteworthy feature observed under salt stress was the excretion of NaCl crystals through salt glands on both the adaxial and abaxial surfaces of the leaves, beginning 15 days after the application of salt stress (Fig. 1B). This salt excretion likely plays a crucial role in osmotic regulation and salt tolerance by minimizing internal sodium toxicity.

Fig. 1

Salinity-induced morphological and anatomical changes in A. littoralis. A Changes in the shoot and root systems of A. littoralis under different NaCl concentrations. B Excretion of salt crystals through the salt glands of A. littoralis. C Root cross-section anatomy of A. littoralis under control and 300 mM NaCl conditions (D). The cell layers are the epidermis (ep), endodermis (en), cortex (co), exodermis (ex), sclerenchyma (sc), the stele (st), central metaxylem (cmx), metaxylems (mx), and aerenchyma (Ae). Leaf cross-section autofluorescence of A. littoralis under control (E) and 300 mM NaCl conditions (F). The structures are as follows: salt gland (Sg), xylem (Xyl), phloem (Ph), sclerenchyma cell (Sc), mesophyll cell (mc), bundle sheath cell (bs), stroma (S), and epidermis cell (ec). Scale bars: 100 μm

Anatomical changes were noticed in the transverse sections of A. littoralis roots and leaves under 300 mM NaCl when compared with those under normal conditions by using the natural autofluorescence of the plant cell walls via epifluorescence microscopy. As shown in Fig. 1C and D, root anatomical changes included enlarged and increased cortical cell layers with the formation of aerenchyma (Table S1). These features may enhance internal oxygen availability and facilitate water and ions storage under saline conditions. Additionally, endodermal cell walls thickened by up to 50% and xylem vessel diameter expanded by approximately 30%, potentially improving nutrient selectivity and water efficiency, respectively (Fig. 1C‒D). Likewise, leaf anatomy was altered under salt stress conditions (Fig. 1E‒F). In the Kranz-type mesophyll, characteristic of C4 plants, cells became more compact, with overall thickness reduced by 23% compared to control plants (Table S1). This compaction likely reduces apoplastic Na⁺ diffusion by limiting intercellular spaces, while preserving photosynthetic efficiency through protected bundle sheath cells. Vascular tissues exhibited marked enlargement, with a 35% increase in xylem vessel diameter and a 70% increase in sclerenchyma layer thickness. These changes are likely adaptive, enhancing water transport and supporting osmotic balance through increased mechanical support and water storage capacity (Fig. 1E‒F). Taken together, the observed morphological and anatomical adaptations in A. littoralis under salt stress suggest a complex, multifaceted strategy for coping with high salinity. These include growth modulation, salt exclusion, structural reinforcement, and physiological regulation to maintain function under adverse conditions.

Chromosome content of A. littoralis

To investigate the genome structure of A. littoralis, the chromosome number was surveyed via DAPI counterstaining. The species was found to possess 10 chromosomes, with an average chromosome size of 1.4 μm ± 0.2 (Fig. 2A). In comparison, Oryza sativa (2n = 2x =24) and Triticum turgidum (2n = 4x = 28) exhibit higher chromosome numbers and larger chromosome sizes (Fig. 2B and C). This indicates that A. littoralis has a relatively compact genome characterized by fewer and smaller chromosomes.

Fig. 2

Chromosome karyotype analysis of Aeluropus littoralis, Oryza sativa, and Triticum turgidum

Genome sequencing and assembly

The raw sequencing reads generated from the A. littoralis genome using the BGISeq-500 platform were processed to remove adapter sequences and low-quality reads. After quality filtering, a total of 27.48 Gb of high-quality paired-end short reads were collected (Table S2). This dataset represents approximately 78.5× coverage of the A. littoralis genome, based on an estimated genome size of ~349 Mb [22]. Additionally, whole-genome sequencing of A. littoralis was performed on the PacBio Sequel II platform in HiFi CCS (circular consensus sequencing) mode. After adapter removal and data filtering via SMRT Link v10.1 software, 518.56 Gb of clean subreads remained, comprising 35,051,483 subreads with a mean subread length of 14.79 kb and a maximum length of 432.8 kb. From these, 1,069,261 HiFi reads were obtained, totaling 18.15 Gb, which corresponds to approximately 51.86× genome coverage. The HiFi reads had a mean length of 17.6 kb and a maximum read length of 49.2 kb (Table S3).

The genome size of A. littoralis was estimated via k-mer-based analysis with the following formula: G = k-mer number/average k-mer depth. The 21-mer frequency distribution analysis revealed that the estimated genome size was ~360 Mb, with heterozygosity rates of approximately 1.88% and a 37.89% repeat content (Fig. S2). These results are consistent with previous estimated genome size ~349 Mb based on flow cytometry [22, 44].

The Aeluropus littoralis genome was assembled and assessed via the QUAST quality assessment tool (version 5.0.2). The final assembly spanned 348 Mb across 4,078 contigs, with sizes ranging from 23.044 kb to 5.104 Mb (Table 1). Based on k-mer analysis and flow cytometry, the assembly covers approximately 96% and ~99% of the estimated genome sizes of ~360 Mb and ~349 Mb, respectively. The assembly metrics included a contig N50 of 178.06 kb and an overall GC content of 45.58%. The genome completeness, as evaluated by BUSCO, reached 88.2%, with 85.3% representing single-copy orthologs and 2.9% duplicated (Fig. 3A). These findings indicated a high level of completeness.

Table 1 Summary statistics for the assembled Aeluropus littoralis genome

Fig. 3

A. littoralis genome assembly assessment and comparative synteny analysis. A Benchmarking universal single-copy orthologs (BUSCO) evaluation results for the assembled A. littoralis genome, indicating high completeness based on conserved gene content. B Syntenic genomic blocks between A. littoralis and selected monocot species. (Eleusine coracana, Zea mays, Setaria italica, Sorghum bicolor, Oryza sativa, and Triticum turgidum). Blue lines denote indicate conserved syntenic genomic regions

Repeat annotation

Within the hybrid assembled A. littoralis genome, we identified a total of 90.6 Mb of repetitive sequences, constituting 26.03% of the assembled genome (Table 2). Transposable elements (TEs) constitute a significant class of repeats, spanning approximately 39.96 Mb, or 11.48% of the genome (Table 2). Additionally, 46.12 Mb of unclassified repeats were identified, along with smaller proportions of other repeat classes, including satellite DNA (34.81 kb), simple repeats (2.89 Mb), and low-complexity repeats (452.48 kb). LTR elements accounted for 6.71% (23.35 Mb) of the total TE content. Among these, Gypsy and Copia elements were the most prevalent, representing 6.59% of the genome, with 13,888 Gypsy and 6,777 Copia elements identified (Table 2).

Table 2 Classification of repeat elements in A. littoralis genome

Gene prediction and functional annotation of the A. littoralis genome

A total of 35,287 genes were predicted and annotated from the assembled A. littoralis draft genome. The average lengths of the predicted genes and coding sequences (CDS) were 13,560 bp and 954 bp, respectively. Functional annotation revealed that 29,049 genes (82.65% of the total predicted genes) showed significant similarity to sequences in existing databases. Comparative genomic analysis highlighted sequence homology between A. littoralis and several related species, including Eleusine coracana, Setaria italica, Oryza sativa, Zea mays, and Sorghum bicolor. Notably, E. coracana, Z. mays, S. italica, and S. bicolor, all C4 plants, exhibited the highest sequence similarity, with identity levels reaching 96.8%, 87.0%, 86.1%, and 83.07%, respectively (Fig. 3B). The functional classification of predicted A. littoralis genes based on Gene Ontology (GO) terms across the three main categories, biological process (BP), molecular function (MF), and cellular component (CC), revealed several highly enriched terms. In the BP category, key enriched processes included regulation of gene expression (GO:0010467) with 1,769 genes, response to stimulus (GO:0050896) with 1,417 genes, RNA metabolic process (GO:0016070), and nitrogen compound metabolic process (GO:0006807), together encompassing a total of 3,904 genes. In the MF category, highly represented functions included nucleotide binding (GO:0000166; 2,559 genes), metal ion binding (GO:0046872; 1,621 genes), ATP binding (1,328 genes), hydrolase activity (GO:0003676; 2,746 genes), and catalytic activity (GO:0003824; 8,031 genes). For the CC category, the most prevalent components were the membrane (GO:0016020; 3,479 genes), organelle (GO:0043226; 3,363 genes), cytoplasm (GO:0005737; 2,821 genes), and nucleus (GO:0005634; 1,353 genes) (Fig. S3). The classification of the A. littoralis predicted genes on the basis of their involvement in different KEGG (Kyoto encyclopedia of genes and genomes) metabolic pathways revealed that a total of 9,297 genes were mapped across 258 different metabolic pathways. Furthermore, pathways were classified into seven different modules on the basis of the canonical classes of pathway maps in the KEGG database. The most enriched pathways were “metabolism”, “DNA repair”, and “transport of small molecules”, followed by “signal transduction”, “gene expression”, “cell cycle” and “cellular responses to stimulus” (Fig. 4).

Fig. 4

Bar graph of enriched KEGG pathways identified in the genome of A. littoralis

Genome-wide analysis of the A. littoralis A20/AN1 stress associated protein family

The A20/AN1 domain-containing stress-associated protein (SAP) family plays crucial roles in plant growth, development, and responses to environmental stresses [46,47,48]. A genome-wide analysis of the A. littoralis genome led to the identification of 13 AlSAP genes, designated as AlSAP1‒13. Detailed analyses were conducted to examine the exon–intron structures, physicochemical properties of the encoded AlSAP proteins, and their predicted subcellular localizations. Notably, the majority of the AlSAP genes were found to be intronless, with the exception of AlSAP13, which contains a single intron within its coding sequence (Fig. 5A). The domain structure and motif arrangement analyses revealed that eight AlSAP proteins contain both A20 and AN1 domains. In contrast, four other AlSAP proteins possess two AN1 domains but lack an A20 domain, except for AlSAP10, which has a single AN1 domain and no A20 domain (Table 3). Notably, AlSAP12 harbors one AN1 domain along with two A20 domains (Fig. 5B). In silico tertiary structure predictions of A. littoralis SAP family proteins generated using the AlphaFold 3 server revealed that all AlSAP proteins comprise α-helices, β-strands, and random coil elements, with a high degree of structural homology (Fig. 5C). Among these, random coils were the predominant structural component, followed by α-helices and then β-strands. The subcellular localization predictions indicated that most AlSAP proteins are likely localized to the nucleus, except for AlSAP2, AlSAP4, and AlSAP8, which were predicted to be cytoplasmic (Fig. 5D).

Fig. 5

Computational analysis of A. littoralis SAP family members. Gene structures of the A. littoralis SAP genes. A The CDS, upstream/downstream UTR, and introns are illustrated by blue, orange, and black lines, respectively. B The conserved domains of the AlSAP proteins analyzed via the MEME v5.4.1 tool. C Tertiary structure analysis of A. littoralis SAP family members predicted via alphafold 3 server (https://alphafoldserver.com/welcome). D Heatmap of the predicted subcellular localization of A. littoralis SAP family members

Table 3 Information of A20/AN1 zinc finger domain containing proteins in A. littoralis

Phylogenetic analysis of AlSAP amino acid sequences and their orthologs in Eleusine coracana, Setaria italica, Sorghum bicolor, Oryza sativa, Zea mays, Hordeum vulgare, and Arabidopsis thaliana revealed that the amino acid sequences of the 96 SAPs were divided into five distinct groups (G1‒G5) (Fig. 6). The AlSAP members were distributed mainly into groups G2 and G5. However, AlSAP13, AlSAP6, and AlSAP3, which harbored two AN1 domains, were clustered in group G1. Additionally, all the identified groups included SAPs from both monocotyledons and dicotyledons (A. thaliana), except for group G3, which included only SAP members from monocot species, suggesting that SAP proteins are ubiquitous and appeared before species divergence. Moreover, several AlSAP proteins closely clustered with E. coracana SAP proteins, such as AlSAP1/EcSAP1, AlSAP5/EcSAP5 AlSAP6/EcSAP17, AlSAP9/EcSAP9, AlSAP10/EcSAP15, AlSAP11/EcSAP11, and AlSAP13/EcSAP16, which suggests high homology and a close relationship between A. littoralis and E. coracana.

Fig. 6

Phylogenetic relationships of stress-associated proteins (SAPs) from Aeluropus littoralis (Al), Eleusine coracana (Ec), Setaria italic (Si), Sorghum bicolor (Sb), Oryza sativa(Os), Zea mays(Zm), Hordeum vulgare (Hv), Arabidopsis thaliana (At), Cucumis sativus (Cs), Jatropha curcas (Jc), Medicago truncatula (Mt), and Solanum tuberosum (St). A total of 153 SAPs were divided into five groups. The accession numbers of the protein sequences used to conduct this phylogenetic analysis are listed in Table S6

Differential expression analysis of stress-responsive genes

The expression profiles of four randomly selected AlSAP genes under salt (300 mM NaCl) and osmotic (10% PEG 6000) stresses were validated via qPCR technology. The results revealed that all the selected AlSAP genes were regulated under salinity and osmotic stress conditions in the roots, leaves or both tissues of A. littoralis (Fig. 7). Interestingly, the expression level of the AlSAP4 gene was the highest, reaching 110-fold greater than that under control conditions, especially in the root tissue subjected to 24 h of salt stress. In addition, the expression levels of the AlSAP5 and AlSAP11 genes were greater in leaves than in root tissues, and their expression levels reached 4.8-fold and 6-fold greater than those under control conditions following 6 h of osmotic stress treatment (Fig. 7). However, the expression of the AlSAP4 and AlSAP12 genes was upregulated in the roots but slightly changed in the leaves under both salt and osmotic stress. Although the evaluated genes belong to the A20/AN1 domain-containing protein family, they presented differential spatio-temporal expression, which suggests their differential functions and modes of action (Fig. 7)

Fig. 7

Expression profiles of selected AlSAP genes in the leaves and roots of A. littoralis subjected to salinity and osmotic stress. Transcript levels of AlSAP genes were quantified in leaf and root tissues following treatment with 300 mM NaCl (salt stress) and 10% PEG-6000 (osmotic stress). Data represent mean expression values ± SE (n = 3 replicates). Different lowercase letters indicate significant differences (p< 0.05) according to Duncan’s test)

In silico analysis of regulatory cis-elements in the promoters of AlSAP genes

To decrypt the cis-elements contributing to the differential spatio-temporal expression of the analyzed AlSAP genes, the 2,000 bp nucleotide sequence upstream of each selected AlSAP gene start codon was retrieved and analyzed in the PlantCARE database. The detailed information on the cis-acting elements can be found in Table S4. As shown in Fig. 8, the promoter regions of all the selected AlSAP genes contained promoter core sequences, such as TATA-box, CAAT-box and GATA-box. Additionally, the most frequently occurring cis-acting elements include light-, hormone-, biotic and abiotic stress-, and development- responsive elements. The G-box, Sp1, MRE, and GT1 motifs were identified as light-responsive elements, whereas the TGACG-motif, CGTCA-motif, TGA-element, and TCA-element were identified as MeJA-, auxin-, and salicylic acid (SA)-responsive elements, respectively. The ABRE, P-box, and GARE-motif were characterized as ABA- and gibberellin responsive elements, respectively. Promoter analysis revealed several cis-elements related to environmental stress responses, including STRE, DRE, WRE, LTR, and TC-rich elements. The enrichment of the AlSAP4 and AlSAP11 promoters on cis-element related to environmental stress responses support the qPCR results mentioned above and confirm the usefulness of the A. littoralis assembled genome (Fig. 8).

Fig. 8

Analysis of the cis-acting elements in the promoter region of selected AlSAP genes. The 2000 bp upstream start codon of the selected AlSAP genes were analyzed with the PlantCARE database (https://bioinformatics.psb.ugent.be/webtools/plantcare/html/). The cis-acting elements were classified into four major classes: hormone responsive cis-elements, development-related cis-elements, light responsive cis-elements, and stress-responsive cis-elements

Discussion

Halophyte plants displayed remarkable potential source of genes for environmental stress tolerance, especially towards salinity [2, 18, 21]. This tolerance involves several molecular, biochemical, and physiological mechanisms [17, 18]. The majority of molecular salt tolerance mechanisms have been characterized in model glycophytes such as Oryza sativa and Arabidopsis thaliana [49, 50]. The availability of high-quality genome sequences and annotations in halophyte plants will unlock their molecular studies in terms of salt tolerance, improving our understanding and providing valuable insights into their tolerance mechanisms [51, 52]. In the present study, the morphological characterization of A. littoralis changes under salt stress conditions revealed marked modulation of root system architecture associated with changes in aerial parts changes characterized by decreases in internodes length, leaf size, and NaCl crystal secretion. These findings suggest that these morphological changes probably help A. littoralis to cope with excess salt, regulate osmotic balance, minimize the harmful effects of high salinity, and thereby thrive in saline environments. Furthermore, these morphological changes were supported by several anatomical changes in A. littoralis roots and leaves in response to salt stress. In roots, these anatomical changes are marked by an increase in the size of cortex cells, xylem vessels, endoderm thickness, and the formation of aerenchyma in root tissue, which probably improve the availability of oxygen, the ability to store water and ions, and efficient water transport. In leaves, the Kranz anatomy with compact mesophyll cells, enlargement of vascular tissues, and the expansion of sclerenchyma cell layers are markedly revealed anatomical changes, which probably increase the ability of leaves to perform efficient C4 photosynthesis, store water and manage osmotic pressure. Usually, the halophyte species displayed two approaches to avoid the toxic effects of Na⁺ accumulation, including increasing leaf succulence to dilute the Na⁺ concentration and excreting salts via bladders and salt glands [17, 18, 53]. In addition to such approaches, several halophyte species exhibited numerous xeromorphic features (thick cuticles, small leaves, and abundant trichomes), which have contributed to their ability to grow commonly in saline land in many locations worldwide [17, 18]. The morpho-physiological plant response to salinity is accompanied by anatomical changes in the structures of leaves, stems, and roots [54]. Several anatomical modifications in response to salt stress have been previously reported in different plant species, especially grasses [55,56,57]. The anatomical changes in A. littoralis subjected to salt treatments are consistent with those previously reported by Barzegargolchini et al. [58], who revealed the expansion of sclerenchyma tissue around vascular bundles under salinity, which helps A. littoralis to enhance water retention and ion sequestration. Similarly, thick layers of sclerenchyma cells on both sides of leaf vascular bundles were revealed in Aeluropus lagopoides subjected to high salinity [59]. Thus, well-developed sclerenchyma can increase the efficiency of vascular tissues in the presence of salt stress, helping manage nutrient uptake and distribution [54]. Additionally, the development of gas-filled spaces in the root cortex of A. littoralis, called aerenchyma, represents a particularly beneficial type of adaptation that facilitates gas exchange and improves oxygen delivery to A. littoralis roots under saline conditions where oxygen availability may be limited. High levels of aerenchyma in the roots of Imperata cylindrica and Sporobolus arabicus grass plants promoted their root performance under stress conditions and helped them maintain growth and yield despite the adverse effects of salinity [60, 61]. The morphological and anatomical adaptations provide a significant advantage to A. littoralis in saline environments, enhancing its ability to tolerate salt stress. To explore the genomic basis behind these anatomical and morphological adaptations, chromosome content, genome size estimation, assembly and annotation of the nuclear genome were conducted.

A survey of the karyotype features of the Tunisian ecotype A. littoralis revealed that it possesses small chromosomes (2n = 2X = 20) with a base chromosome number of X = 10 and a genome size estimated to be 360 Mb. These findings are in agreement with those of a previous study on Aeluropus species (A. macrostachyus, A. littoralis, and A. lagopoides) reported by Maryam et al. [62], which reported the existence of only one ploidy level (2n = 20) for all the Aeluropus analyzed species. Generally, the four basic chromosome numbers 7, 9, 10, and 12 are the most common basic chromosome numbers in Poaceae [63]. On the basis of the consistency of our results and those of several previously published studies [22, 44, 62], the number of basic chromosomes within the A. littoralis species can be considered stable with x = 10 chromosomes. Owing to its reduced chromosome number with small size, A. littoralis is likely an excellent halophyte model organism for genetic research.

To successfully assemble a high-quality A. littoralis genome, a hybrid approach was adopted using short and long reads from the BGISeq-500 and PacBio HiFi reads, respectively. The assessment of the assembled A. littoralis genome revealed interesting features, including a reasonable contig number (4078), a longest contig length of 5.1 Mb, an acceptable N50 value of 133.77 kb, and a high completeness percentage value of 91.1% on the basis of the presence of orthologous Poaceae genes. Compared with the previously described A. littoralis genome reported by Hashemi-Petroudi et al. [44], our assembled A. littoralis genome was characterized by a reduced contig number (4078/182,747 contigs), and greater predicted gene number (35,147 vs. 15,916 genes) and N50 value (133.77 kb vs. 3.6 kb), which confirms that our A. littoralis assembled genome is of greater quality and that longer contigs contribute to a more complete representation of the A. littoralis genome. Since Hashemi-Petroudi et al. [44] used only short reads collected from the Illumina HiSeq 2500 platform, our hybrid approach, especially the use of PacBio Hifi long reads, contributed to considerable improvements in all the metrics related to the quality of the assembled A. littoralis genome. A total of 35,147 genes were predicted, which is quite similar to the number of genes predicted in several Poaceae species, such as S. italica (34,584 genes), O. sativa (39,387 genes), E. coracana (48,836 genes), Z. mays (49,897 genes), S. bicolor (34,118 genes), and P. halii (33,805 genes).

The synteny inferred on the basis of the predicted genes from the species O. sativa, E. coracana, S. italica, Z. mays, S. bicolor, and T. turgidum highlighted the highest sequence similarity of 96.8% between the A. littoralis and E. coracana sequences, followed by 87%, 86%, 83.07%, 79.19% and 64.18% with the Z. mays, S. italic, S. bicolor, O. sativa, and T. turgidum sequences, respectively. These findings were consistent with previously published phylogenetic analyses based on chloroplast genome sequences, which revealed that the Eleusine coracana and Eleusine indica species from the Euleusininae subtribe are the closest sisters to A. littoralis, which belongs to the Aeluropodinae subtribe [25].

The annotated repeat elements represented 26.03% of the A. littoralis assembled genome and were distributed in different classes of repetitive elements, including retroelements and DNA transposons. Repetitive sequences are considered significant components of plant genomes, and their presence can play a role in genome evolution, gene regulation, and adaptation to environmental stresses [64, 65]. The fraction of repeated elements in the assembled genome of A. littoralis is quite comparable to that annotated in monocot plants with small genomes, including O. sativa (39%) [66] and Oryza coarctata (25.5%) [67]. However, given the positive correlation between plant genome size and the number of repeat elements in the genome, it is more appropriate to compare repeat element features across plant genomes. For instance, Gypsy LTR elements are considered major determinants of genome size expansion in plants [68]. Thus, repeat element analysis revealed that the Gypsy LTR elements (5.05%) are more represented than the Copia LTR elements (1.54%), which is consistent with most grass genomes, including O. sativa [66], Z. mays [67], S. bicolor [69], and S. italic [70]. Future studies are recommended to explore the contribution of the annotated repeat elements to the gene regulation and adaptation of A. littoralis to environmental stresses.

On the basis of homology searches against the Swiss-Prot, InterProt, and NCBI non-redundant protein databases, 82.65% (29,049 genes) of the total predicted genes were successfully annotated. The functional categorization of the predicted genes via gene ontology (GO) analysis revealed that the enriched GO terms involved several biological processes, including the regulation of gene expression, response to stimulus, and nitrogen metabolic processes. The markedly enriched molecular functions in which the predicted genes are involved were nucleotide binding, metal ion binding, ATP binding, hydrolase activity, and catalytic activity. Mapping of the A. littoralis predicted genes into different KEGG metabolic pathways revealed that these genes were involved mainly in metabolism, DNA repair, signal transduction, gene expression, and cellular responses to stimulus. Plants are constantly exposed to environmental stresses and to overcome these conditions they resort to induce the expression of several stress gene families [19]. Although the key role of A20/AN1 stress-associated proteins (SAPs) in environmental stress tolerance has been described for many plant species [46,47,48, 71, 72], genome-wide information on this gene family in A. littoralis is still lacking. Recently, the AlSAP gene encoding an A20/AN1 stress-associated protein was isolated from A. littoralis, and its overexpression enhanced the salt and osmotic stress tolerance of transgenic tobacco, rice and durum wheat lines [26, 27, 30,31,32,33, 46]. Through genome-wide analysis of the assembled A. littoralis genome sequence, thirteen AlSAP genes were identified. The genomic features, physiochemical properties and the subcellular localization of the encoded AlSAP proteins were also determined. Notably, while most AlSAP proteins were computationally predicted to localize to the nucleus, three members‒AlSAP2, AlSAP4, and AlSAP8‒exhibited cytoplasmic localization. These findings align with previous reports on SAP proteins in other species. For instance, Wang et al. [73] reported that castor bean SAP-GFP fusions localized to multiple cell compartments, including nucleus, cytoplasm, and cytomembrane, in tobacco epidermal cells. Similarly, Sidra et al. [74] noted the nuclear and cytoplasmic localization for nine almond PdSAP proteins. The number of A20/AN1 gene members varies widely across plant species but generally ranges from 10 to 20 genes [46, 48]. The monocotyledonous species closest to A. littoralis have a quite similar numbers of A20/AN1 gene members, with 10 members in Z. mays [71], 12 members in S. italic [75], 18 members in both rice and S. bicolor [76, 77], and 17 members in H. vulgare [78]. Several studies have reported that A20/AN1 zinc finger gene family members are involved in environmental stress responses, including osmotic and salt stress [46, 48, 75,76,77]. The monitored transcript level accumulation of four selected AlSAP genes (AlSAP4, AlSAP5, AlSAP11, and AlSAP12) revealed their quick and differential regulation in the roots and leaf tissues of A. littoralis subjected to salt and osmotic treatment. The response to environmental stresses requires the presence of cis-elements in the promoter regions of genes to drive their transcription [79]. A survey of the regulatory cis-elements in the promoters of the selected AlSAP genes revealed the presence of light- (G-box, Sp1, MRE, and GT1 motifs), hormone- (TGACG-motif, CGTCA-motif, TGA-element, TCA-element, ABRE, P-box, and GARE-motif), and biotic/abiotic stress (STRE, DRE, WRE, LTR, and TC-rich element)-responsive elements. The large number of cis-elements related to environmental stress responses in the promoter regions of the AlSAP4 and AlSAP11 genes supports their rapid and highly differential expression in response to salinity and osmotic stress. The presence of hormone-responsive cis-elements, including GA-responsive cis-elements (P box), auxin-responsive element (TGA-element), salicylic acid responsive element (TCA-element), and MeJA-responsive element (TGACG- and CGTCA-motif) suggests that the analyzed AlSAP genes are potentially regulated by the above mentioned phytohormones and may play roles in crosstalk between defense signaling pathways and growth signaling pathways. Furthermore, the identification of development-related cis-elements, such as the endosperm-related elements (GCN4 motif) and meristem-related motifs (CCGTCC-motif and CAT box), implies that these analyzed AlSAP genes may exhibit tissue-specific expression patterns during plant growth and differentiation. These findings align with previous functional studies on SAP gene promoters and further supporting the stress-responsive nature of SAP genes regulation [27, 80,81,82]. Moreover, these evidences suggest differential functions and modes of action of the analyzed AlSAP genes and confirm the usefulness of the assembled A. littoralis genome.

In future work, to enhance the quality of our draft genome assembly, we plan to generate Hi-C sequencing data. This will allow us to bridge gaps between contigs and reconstruct the chromosome-scale architecture of the A. littoralis genome. Hi-C technology leverages three-dimensional genomic organization to identify centromeric and telomeric regions, characterize structural breakpoints, and delineate individual chromosomes within an assembly [83]. For example, Hi-C data was used in assembling the cucumber genome into its seven expected pseudo-chromosomes [84]. The genomic synteny analysis revealed high sequence similarity between predicted genes from the A. littoralis draft genome and those of Z. mays, S. bicolor, O. sativa, and T. turgidum, with identity levels of 87%, 83.07%, 79.19%, and 64.18%, respectively. These findings highlight the need for further investigation using RNA sequencing (RNA-seq) of A. littoralis roots and leaves, which will help elucidate gene expression differences and identify key salt stress-responsive genes and their regulatory elements. Such investigation is a crucial first step toward identifying orthologous genes in major cereal crops such as rice, maize, wheat, and sorghum. This knowledge could pave the way for utilizing transgenic approaches or CRISPR/Cas9 genome editing to enhance the tolerance of these cereals to adverse environmental conditions. Similar integrated genomic and transcriptomic analyses have been used to clarify the genomic basis of salt tolerance in zoysiagrass (Zoysia spp.), a warm-season turfgrass known for its tolerance to salt, drought, and heat stress [85].

Conclusions

This study presents valuable insights into the adaptation mechanisms of A. littoralis, a halophyte grass, to saline environments. It highlights both morphological/anatomical adaptations and offers the first high-quality nuclear genome assembly for this species. This assembled genome serves as a robust platform to explore the genetic basis of salinity tolerance in halophytes, which could identify critical genes for improving salinity and drought tolerance in cereals. Additionally, the assembled A. littoralis nuclear genome combined with the previously published chloroplast genome provides a broader understanding of the plant’s evolutionary history within the Poaceae family, particularly the C4 plants. The use of PacBio HiFi and Hi-C sequencing technologies is suggested for future research to generate a chromosomal-level genome.

Methods

Plant materials, salt stress response assays, and histological analyses

Seeds of an Aeluropus littoralis ecotype were collected from salt marshes near “Sfax-Tunisia” (34°37'57.0"N 10°38'47.0"E) and cultivated for several years in a greenhouse for self-pollination at the Centre of Biotechnology at Sfax (Tunisia) and the College of Food and Agricultural Sciences (King Saud University-Saudi Arabia). The specimen used in this study was deposited under voucher number 51,203 in the herbarium of the college of food and agriculture sciences, King Saud University. Taxonomic identification was confirmed by Prof. Dr. Abdulaziz Assaeed, who is affiliated with the College of Food and Agriculture Sciences, King Saud University. Seeds from this ecotype were surface sterilized with 70% ethanol and substantially germinated in MS media. Seven-day-old A. littoralis seedlings were transplanted into MS media and MS media supplemented with 150 mM, 300 mM, 450 mM, or 550 mM NaCl and grown hydroponically for 120 days under controlled conditions at 25 °C with a 16-hour photoperiod. The morphological changes were documented with a Canon EOS 750D camera.

At the end of the salt stress response assays, the root region located 3 cm from the root apex and the youngest fully expanded leaf of A. littoralis plants cultivated in MS and exposed to 300 mM NaCl were removed and fixed with paraformaldehyde (PFA) fixative solution. The fixed A. littoralis roots and leaves were embedded in a 2.5% agar block (Euromedex, Ref. LE-8200-B). Root and leaf cross sections (50 µm in thickness) were obtained with a Microm HM 650 V microtome (Thermo Fisher Scientific, Germany) via the following parameters: frequency, 50 Hz; amplitude, 0.9; and speed, 30, as previously described by Lartaud et al. [86]. Afterward, the autofluorescence of the cell walls of transverse sections was observed under a Leica DM 6000 epifluorescence microscope (Leica Microsystems, Germany), and images of the best sections of roots and leaves were acquired with a color Retiga 2000R camera (Qimaging, Canada).

Chromosome preparation from A. littoralis, O. sativa, and T. durum

Roots were harvested from plants cultivated in a hydroponic system. The samples were treated with 0.04% hydroxyquinoline for 4 h, fixed for 48 h at 3:1 (ethanol:acetic acid) and stored in 75% ethanol at 4 °C. The fixed roots were rinsed twice in H₂O for 10 min each, treated in 0.25 N HCl for 10 min, rinsed for 10 min in H₂O and placed in digestion buffer (0.01 M citrate buffer pH 4.5, 0.075 M KCl) for 10 min. The root tips were cut and placed in an enzyme mixture (1% Onozuka R-10 cellulase, 1% Y-23 pectolyase and 1% cythoelicase in digestion buffer) in a micro tube at 37°C for approximately 3 h (the time varied with the size of the roots). The root tips were then rinsed in H₂O and spread on a slide with a 3:1 ratio of ethanol: acetic acid. The quality of the slides was controlled by microscopic observation via phase contrast. The slides were mounted in Vectashield antifade solution with DAPI (Vector Lab). The slides were examined with a Leica DMRAX2 fluorescence microscope with a triple band filter for simultaneous visualization of color, i.e., DAPI (blue). The images were acquired with a cooled high-resolution black and white CCD camera.

DNA extraction, DNA library preparation, and library sequencing

High-molecular-weight genomic DNA (HMW-gDNA) was extracted according to the Murray and Thompson CTAB method [87] from 2-month-old A. littoralis plants grown hydroponically [32] under greenhouse conditions as previously detailed by Ben Romdhane et al. [25]. The amount and quality of the extracted genomic DNA were assessed via spectrophotometric and electrophoresis methods, respectively [25]. Following a clean-up step using AMPure® PB beads (Beckman Coulter) followed by elution with elution buffer, the quality of the genomic DNA sample was inspected using an Agilent 2100 Bioanalyzer (Agilent Technologies). The A. littoralis genomic DNA sample was subsequently divided into two aliquots in order to be subsequently fragmented and used for BGISeq-500 and HiFi sequencing library construction.

An aliquot of the cleaned A. littoralis genomic DNA sample was mechanically fragmented using Qsonica Q800R2 sonicator. The resulting fragments were ligated with BGISEQ-500 adapters as described by Mak et al. [88]. Following the purification of the library using PB binding buffer (Qiagen, Hilden, Germany) and Monarch DNA Cleanup Columns (New England Biolabs, Massachusetts, USA), the library was eluted in 40 μL of buffer EB (Qiagen, Hilden, Germany). Post-PCR amplification, as an extra purification step, was performed via the AMpure XP system (Agentcourt, Beckman Counter, Indianapolis, USA) with a 1.8 × bead: library ratio to remove any persisting primer dimers or other molecules with a fragment size of <100 bp. The concentration and fragment size of the library were quantified using a high-sensitivity DNA assay kit with the 2100 Bioanalyzer system (Agilent Technologies, Waldbronn, Germany). The 150-bp PE sequencing of the libraries was carried-out using the BGISeq-500 sequencing platform of the Beijing Genomics Institute (BGI-Shenzhen, Wuhan, China). The second aliquot of the cleaned A. littoralis genomic DNA sample was subjected to fragmentation and purification steps as recommended by PacBio and previously detailed by Ben Romdhane et al. [25]. Two libraries were constructed, and A. littoralis whole-genome sequencing (WGS) was subsequently conducted via the DNA Link Sequencing Lab (DNA Link Inc, Seoul, Republic of Korea) using the PacBio Sequel II platform in HiFi circular consensus sequencing (CCS) mode.

A. littoralis hybrid de novo genome assembly and assessment

The collected BGISeq-500 raw reads were processed to remove adapter sequences and low-quality sequences using Cutadapt program [89]. The generated PacBio HiFi reads were pre-processed via SMRT Link software v12.0 (PacBio) and subsequently subjected to a taxonomic classification step with a 0.05 confidence rate using the Kraken2 taxonomic classifier [90] with a prebuilt MiniKraken v2 database including bacteria, archaea, and virus genomes RefSeq. Finally, all reads that were classified as bacterial, archaeal, or viral contaminants were eliminated. The clean collected PacBio HiFi reads (≥17 kb, Q ≥ 20) were subsequently assembled with short reads (100 bp, Q≥ 20) into a contig-level assembly using CLC Genomics Workbench V22.0 software with following parameters: Error correction= iterative, minimum overlap length= 3000 bp, bubble size= auto-detected, mismatch cost= 2, and insertion/deletion cost= 3. After the initial assembly was complete, the haplotigs purging step was performed via purge pipeline v1.1.1 to improve genome continuity by reducing heterozygous duplication [91]. The assembled whole A. littoralis genome was assessed for contiguity and completeness via the Quast tool [92] and the benchmarking universal single-copy orthologs (BUSCO) [93] completeness score, respectively. The aforementioned analysis employed BUSCO v5.2.2 and utilized the poales odb10 dataset.

The A.littoralis genome size was estimated via the k-mer method base on sequencing data generated via PacBio Sequel II in HiFi mode. The clean reads were subjected to 21-mer frequency distribution analysis using a Jellyfish k-mer counter [94] and the GenomeScope program [95].

Repeat annotation

A de novo repeat library based on the assembled A. littoralis genome was generated using RepeatModeler v2.0.1. This library was subsequently merged with the green plant repeat library from the Repbase database version 29.04 (https://www.girinst.org/repbase/) and used to detect repetitive elements via RepeatMasker v4.1.5, and subsequently masked the A. littoralis assembled genome.

Gene prediction and functional annotation

Protein-coding gene prediction was conducted using a combination of homology-based and ab initio prediction via GeMoMa [96] and Augustus v3.5.0, respectively. The gene models of Setaria italica (GCA_000263155), Aegilops tauschii (GCA_002575655), Oryza sativa Japonica Group (GCA_001433935) and Sorghum bicolor (GCA_000003195.3) were downloaded from the European Nucleotide Archive (https://www.ebi.ac.uk/ena/browser/home).

Blast of the protein-coding genes in the Swiss-Prot, eggNOG, and the NCBI non-redundant protein (NR) databases allowed functional annotation. InterProScan v5.35 was used to annotate protein domains and motifs, and the matching InterPro entries were used to assign gene ontology (GO) terms.

Analysis of the AN1/A20 domain-containing gene family

Genome-wide analysis of AN1/A20 zinc finger domain containing proteins

The stress associated proteins (SAPs) containing AN1/A20 zinc finger domains from Oryza sativa, Aegilops tauschii, Sorghum bicolor, and Setaria italica were retrieved from the Ensemblplants database (https://plants.ensembl.org/index.html) and subsequently used as queries to identify A. littoralis homologs via the tBLASTn tool. The candidate A. littoralis SAP proteins were then screened using SMART (http://smart.embl-heidelberg.de/) to reveal the existence of the PF01754 zinc finger A20 domain and PF01428 zinc finger AN1 domain. Then, the theoretical isoelectric point (pI) and molecular weight (MW) of the all discovered AlSAP proteins were subsequently determined via ExPASy tools (https://www.expasy.org/). The A20/AN1 motifs conservation was established using the MEME v5.4.1 tool (https://meme-suite.org) [97]. Protein motifs and gene structure were visualized via TBtools software [98].

The WoLF PSORT (https://www.genscript.com/wolf-psort.html) and CELLO v.2.5 (https://cello.life.nctu.edu.tw/) were utilized to estimate the subcellular location of the identified AlSAP proteins based on their functional motifs and amino acid composition. Additionally, the three-dimensional structures of AlSAP proteins were modeled by the AlphaFold 3 server (https://alphafoldserver.com).

The AN1/A20 zinc finger protein sequences of Oryza sativa [99], Setaria italica [75], Arabidopsis thaliana [100], A. littoralis, Sorghum bicolor, Zea mays, Eulecine coracana and Hordeum vulgare [101], Cucumis sativus, Jatropha curcas, Medicago truncatula, and Solanum tuberosum were retrieved to develop a phylogenetic tree using maximum likelihood method with 1000 bootstrap replicates, as implemented by MEGA12 [102].

To determine the putative cis-acting elements into four selected AlSAP gene promoters (AlSAP4, AlSAP5, AlSAP11, and AlSAP12) and their possible involvement in stress regulation and plant development, the 2,000 bp upstream of the transcription start site as promoter regions were used as input to perform the analysis via the PlantCARE database (http://bioinformatics.psb.ugent.be/webtools/plantcare/html/).

Transcript accumulation analysis of AN1/A20 zinc finger proteins

Total RNA was collected from frozen A. littoralis leaf and root samples via an E.Z.N.A Kit (Omega Bio-Tek). For RT-qPCR, one μg of total RNA per sample was employed to synthesize first-strand cDNA using a QuantiTect RT-Kit (Qiagen). The transcript accumulation of the AN1/A20 domain containing genes was monitored via RT-qPCR reactions as previously outlined by Ben Romdhane et al. [33]. Relative expression was determined by the comparative threshold cycle (2^−ΔΔCT) method [103], using the 26S rRNA-based internal control gene (Table S6).