NOD-like receptor repertoire in the chromosome-level genome of the demosponge Dysidea avara (Schmidt, 1862)

Porifera, one of the earliest diverging metazoans, have shown a surprisingly complex immune repertoire. However, most information to date is based on de novo transcriptome assemblies, limiting our knowledge regarding the presence and evolution of poriferan immune repertoire. Here, we generated the chromosome-level genome of the demosponge Dysidea avara, a target species in studies on symbiosis and differential expression of immune genes. We examined the presence and the number of common immune protein domains in the annotated genome of D. avara, and we further focused on NOD-like Receptors (NLRs), which are one of the most expanded immune receptors in Porifera according to previous reports on draft genomes and transcriptome assemblies. Dysidea avara has a 575 Mb genome with N50 41Mb, 162 scaffolds, and 15 chromosomes. We additionally recovered 37 sequences corresponding to microbial genomes, including complete bacterial and viral genomes. Based on the presence of conserved domains, we detected a large number of immune receptors and other immune genes in D. avara genome, such as 14 TIR, 39 CARD, 128 DEATH, and 230 NACHT domain-containing genes. Based on their architecture, we identified a large expansion of bona fide NLRs (i.e., 126 NACHT+LRR domain-containing genes); of which, 20 included a N-terminal CARD domain (NLRC), and 25 included a N-terminal DEATH domain (NLRD). In D. avara, the different NLR categories (i.e., NLRX, NLRC and NLRD) formed distinct phylogenetic clusters, while the NLR phylogenetic analysis across sponge chromosome-level genomes indicated that NLRs were mainly grouped by species rather than category. The NLRX category was the most expanded, while the NLRC category was absent in 7 out of 11 studied sponge genomes. These observations indicate that the diversification of NLRs in sponges, most likely derived from the ancestor NLRX, responds to species-specific selective pressures related to their immunity. This is the first study characterizing sponge NLR diversity in a chromosome-level genome, enhancing our knowledge of NLR evolution in the ancient phylum Porifera.

1 Introduction

Porifera (sponges) is one of the earliest phyla in the evolution of metazoans (1). They are sessile organisms, dating back to 500–700 million years ago (2–4) with more than 9,300 species recorded to date (5). The phylum Porifera is divided into four different classes (Hexactinellida, Calcarea, Homoscleromorpha, and Demospongiae (6–9)), representing a crucial group for studying the evolutionary transition from unicellular to multicellular life and the emergence of metazoan traits (10). More recently, sponges have emerged as a prominent group to understand the evolution of animal-microbe interactions. Like most metazoans, sponges maintain intimate relationships with a large variety of microbes (11). It is noteworthy that sponges preserve a consistent, well-defined microbiome (12–14) while continuously filtering thousands of liters of seawater per day and encountering a wide array of foreign and potentially harmful bacteria (15). The maintenance of a stable microbiome, together with the possible recognition of distinct microbes in sponges is thought to be orchestrated by their immune system (16).

The innate immune system plays an essential role in interkingdom cross-talk and the maintenance of a healthy symbiosis between animals and microbes (17–19). Increasing evidence shows the role of host Pattern-Recognition Receptors (PRRs) in the recognition of beneficial symbionts and host-microbiome communication (16, 20). Several PRR families present recognizable architectures composed of at least one conserved protein domain that is present in homologs across the different animal groups (21). The nucleotide-binding and leucine-rich repeat receptors (NLRs) are among the best characterized family of PRRs (22). Interestingly, NLRs are especially diversified in early diverging metazoans, such as cnidarians and sponges (23–26) and in other marine invertebrates (27), and are hypothesized to play a role in invertebrate immune specificity (28, 29). NLRs are recognized by the conserved NACHT protein domain (30, 31), combined with variable C-terminal LRR domains (30). This architecture defines bona fide NLRs. Their canonical structure, as described in vertebrates, is complemented by an N-terminal with CARD, Death, Pyrin, or BIR domain (30, 31). NLRs in vertebrates recognize microbes via the C-terminal LRR, which bind microbial patterns (e.g., lipopolysaccharide (LPS), peptidoglycan (PGN), peptides) (32, 33). Then, the N-terminal domain participates in protein-protein interactions and is responsible for activating further downstream pathways related to inflammation during pathogen infections (34). NLRs in invertebrates seem to have a similar role, being activated in defense mechanisms (23, 35). Besides, NLRs have also been considered to play a role in differentiating mutualists from pathogens and in symbiosis establishment, for example in the intestinal microbiota in humans (36). Similarly, the expanded NLR repertoire in the coral Acropora digitifera have also been proposed to be devoted to interaction with their obligate dinoflagellate endosymbiont (24).

Previous studies based on draft genomes and transcriptomes have indicated a complex innate immune repertoire in sponges (11, 26, 28, 37–41). Sponges possess a relatively large expansion of bona fide NLR receptors, which might indicate a high level of specificity to the large range of microbes they encounter (28). Functional validation of NLR roles in sponges is limited, but bulk transcriptomics in the sponge Dysidea avara revealed the differential expression of NLR genes in response to a “cocktail” of microbial-associated molecular patterns (MAMPs) (26) and when incubated with sponge-associated vs food bacteria (42). Building on this previous knowledge, our focus here was the species D. avara. Dysidea avara is an encrusting sponge of the class Desmospongiae from the Atlantic-Mediterranean region, found in rocky sublittoral habitats, at a maximum of 80m depth (43, 44). It is characterized as Low Microbial Abundant (LMA) sponge because it harbors a low quantity of symbiotic bacteria, compared to other species (45). This species has received considerable attention due to its potential pharmaceutical value. Dysidea avara produces avarol and derivatives of broad-spectrum activity, including antitumor, antipsoriatic, and antileukemic properties (46, 47). Moreover, the physiology of D. avara and its role in benthic-pelagic nutrient transfer has been extensively studied, connecting the metabolism of the sponge and its microbiota with ecosystem processes (48–50). This species is one of the few sponge species in which larva settlement, juvenile development and metamorphosis have been achieved in vitro (51). Due to all the above, D. avara emerges as a model species for studying sponge-microbe interactions and the evolution of the innate immune system (26, 42, 52). Here, we generated a chromosome-level assembly from this species with PacBio and Hi-C technologies with the aim of identifying the presence/absence and expansion of most common metazoan immune related domains and further exploring the NLR receptor expansion and structure in this species. We further conducted a phylogenetic analysis to understand the NLR evolution and expansion within D. avara and in relation with other sympatric sponge species. In comparison to previous investigations based on de novo transcriptome assemblies (e.g. 26, 42), the use of such a chromosome-level genome assembly in NLR investigation is expected to give accurate information of the full NLR composition (and not only the expressed fraction) and the exact position of those NLR genes on the chromosomes of D. avara. This information can be further used to investigate evolutionary traces of the different NLR families in a broader context and among sponges but also across other animal phyla.

2 Materials and methods

2.1 Specimen collection

An individual of D. avara (Figure 1A) was originally collected in April 2019 at l’Escala, Girona, Spain (42.1145863 N, 3.168486 E) and kept in a flow-through system with direct intake of seawater in aquaria facilities at ICM-CSIC (Barcelona, Spain) for about 1 month. Tissue sample from the sponge was then collected, cleaned from any epiphytes and snap-frozen in liquid nitrogen, followed by storage at -80 °C until further processing. We performed a sponge cell separation by differential centrifugation prior to DNA purification, adapting the protocol by Wehrl et al. (53). This step was meant to enrich the sponge cell fraction of the sample, reducing potential bacterial contamination in DNA extracts. In short, ca. 15 mL of tissue volume was cut into small pieces and incubated in calcium-magnesium-free artificial seawater with EDTA (CMFASW-E) to favor cell dissociation and gentle homogenization, and then filtered through 100 µm Nitex. The resulting suspension was centrifuged at 400x g for 10 min at 4 °C to recover the pellet with the sponge cells. The supernatant was discarded, and the pellet was re-suspended in CMFASW-E and centrifuged again to remove any remnant debris. The resulting pellet was immediately processed for DNA extraction.

2.2 Genome assembly and annotation

2.2.1 High-molecular-weight DNA purification

High-molecular-weight DNA was extracted with the NucleoBond^® HMW DNA kit (Macherey-Nagel Gmbh & Co, Germany) following the manufacturer’s instructions, by adding an extra wash step with solution H4. The quality, quantity, and size of DNA extracts were assessed by NanoDrop 2000c Spectrophotometer (peolab, Germany), Qubit 2.0 (Life Technologies, Carlsband, CA), and by pulse-field gel electrophoresis in 0.75% agarose gel in 1x Lönning buffer with SYBRSafe DNA stain (1:100000; Sigma-Aldrich).

2.2.2 Library construction and draft genome sequencing

Long read sequencing was performed using the PacBio Sequel II platform at CCGA sequencing facility in Kiel. Three consecutive SMRT Cells were run, obtaining a total of 32.1 Gb of CCS reads. The IPA v1.8 (https://github.com/PacificBiosciences/pbipa) assembler was run with local mode (including polishing, purging haplotigs, and phasing) to assemble together the reads from all the three libraries and to generate the D. avara draft assembly. This assembly was used as input for Hi-C sequencing as described below.

2.2.3 Dovetail Omni-C library preparation and sequencing

Frozen material from the same D. avara individual was sent to Cantata company (previous Dovetail Genomics) for generating Omni-C^® libraries and scaffolding the draft genome with the HiRise™ software. For each Dovetail Omni-C library, chromatin was fixed in place with formaldehyde in the nucleus. Fixed chromatin was digested with DNase I and then extracted. Chromatin ends were repaired and ligated to a biotinylated bridge adapter followed by proximity ligation of adapter adapter-containing ends. After proximity ligation, crosslinks were reversed and the DNA purified. Purified DNA was treated to remove biotin that was not internal to ligated fragments. Sequencing libraries were generated using NEBNext Ultra enzymes and Illumina-compatible adapters. Biotin-containing fragments were isolated using streptavidin beads before PCR enrichment of each library. The library was sequenced on an Illumina HiSeqX platform to produce ~ 30x sequence coverage.

2.2.4 Scaffolding with Omni-C HiRise

The input de novo assembly and Dovetail OmniC library reads were used as input data for HiRise, a software pipeline designed specifically for using proximity ligation data to scaffold genome assemblies (54). Dovetail OmniC library sequences were aligned to the draft input assembly using bwa (https://github.com/lh3/bwa). The separations of Dovetail OmniC read pairs mapped within draft scaffolds were analyzed by HiRise to produce a likelihood model for genomic distance between read pairs, and the model was used to identify and break putative misjoins, to score prospective joins, and make joins above a threshold. Genome completeness was assessed with BUSCO v5.3.0 (55), using both the eukaryota_odb10 and the metazoa_odb10 databases to interrogate the assembly for the presence of common single-copy orthologs. A prediction of genome size for D. avara was obtained using GenomeScope software (v2.0) (56) and all generated CCS reads.

2.2.5 Genome annotation

The annotation of the generated chromosome-level genome of Dysidea avara was conducted with the genomeannotator pipeline developed in-house (https://github.com/marchoeppner/genomeannotator). It is an automatic genome annotation pipeline based on ab-initio gene prediction as well as experimental and additional model hints from multiple possible sources, using the Nextflow workflow language (57). We also used a subset of RNA seq data from a previous publication (26), (Experiment ArrayExpress accession E-MTAB-6757; ENA run accessions: ERR2560048; ERR2560049; ERR2560053; ERR2560055; ERR2560056).

2.3 Immune Pfam domain identification

2.3.1 Search of immune Pfam domains in annotated D. avara genome

Specific domains of the metazoan immune repertoire were selected, based on the review by Buckley et al. (21), and searched in the translated genome of D. avara. We used the respective PFAM database entry from the Interproscan page (https://www.ebi.ac.uk/interpro/entry/pfam/) for each domain. The search in the predicted gene models was conducted with hmmerscan v.3.4 (http://hmmer.org/) with default parameters. We further checked those sequences in SMART (http://smart.embl-heidelberg.de/) keeping only the sequences above threshold that included the target Pfam domain as visible or overlapping with other domain(s) within the sequence, as determined by SMART. We expanded our analysis to 10 other translated chromosome-level sponge genomes, recently generated during Aquatic Symbiotic Project by the Wellcome Sanger Institute (ASG), (Supplementary Table 1). We included 4 demosponges and 1 Homoscleromorpha that can all occur in sympatry with D. avara in the Mediterranean Sea. We also included 2 additional demosponges that have been target species of studies on immunity (Halichondria panicea and Amphimedon queenslandica), plus 1 freshwater demosponge (Ephydatia muelleri) and 1 Hexactinellida (Aphrocallistes vastus). The ASG recently generated a chromosome-level genome of our target species, D. avara, which was also included in our analysis for comparison. Finally, we reported the number of genes containing the target Pfam domains for each species.

2.3.2 Search and characterization of D. avara annotated genes coding for NLRs

We extracted the protein sequences of those genes that had similarity with NLRs from the automated eggNOG annotation of our genome. Additionally, we used the NACHT containing proteins retrieved from the previous hmmer search (section 2.3.1) in order to include more divergent D. avara NLRs. We checked manually for their domain architecture in SMART in Genomic mode (http://smart.embl-heidelberg.de/) (58). To illustrate the domain architecture of the candidate NLR proteins, we used the Biorender.com online software. We also checked the positions of those genes on the chromosomes of the D. avara genome, using the RIdeogram v0.0.2 in R (59).

2.3.3 Phylogeny of poriferan NLRs

For determining the phylogeny of D. avara NLR genes in the context of NLR repertoires in other sponge species, we focused on their most conserved domain: the NACHT domain (30, 31). For this analysis we did a different approach in order to extract comparable information for all studied species. We first constructed a custom hmmer profile from the alignment of NACHT domains retrieved in a total of 28 sequences from different phyla of metazoans, placozoans, and other sponges from NCBI (Supplementary Table 2A; Supplementary File 1). We performed alignments with MUSCLE algorithm (60) in Geneious software (61), and selected the part of the alignment that was better aligned among all phyla (Supplementary File 2A). We chose MUSCLE because of its high accuracy on amino acid datasets, without incurring substantial computational cost. Checking in SMART, we further confirmed that this selected alignment was assigned to a NACHT domain. The hmmer profile was used to retrieve NACHT-containing proteins from the translated genome of D. avara (Supplementary Table 2B; Supplementary File 1) and the 10 other sponge translated chromosome-level genomes (Supplementary Table 1; Supplementary Table 2C; Supplementary File 1). The domain architectures of the retrieved NACHT-containing proteins were further checked in SMART, and only those corresponding to bona fide architectures (that is those combined with leucine-rich repeats) were selected for further analysis (Supplementary File 1). In order to generate a phylogenetic tree of bona fide NLRs in D. avara, we combined the results obtained from the hmmer search and the direct search on automated annotation (see section 3.2.3) (Supplementary File 1; Supplementary File 2B). For the NLR phylogeny across sponge species, we also included as outgroup the NLR sequences of the cnidarian species Nematostella vectensis (Supplementary Table 2A; Supplementary File 1; Supplementary File 2C). Both phylogenetic trees were constructed based on the NACHT domain alignment from all the NLR sequences either within D. avara or across all sponge species (Supplementary File 1; Supplementary Files 2B, 2C). The phylogenetic trees were generated by maximum likelihood in RAxML v.8 (62) with the raxmlGUI platform v2.0.1 (63), using GTR Bootstrap model for proteins and an estimated gamma shape parameter. The node support was calculated with a thorough bootstrap algorithm and 100 independent searches. The tree was further processed with figtree v.1.4.4. and the online tool iTOL (https://itol.embl.de/).

3 Results

3.1 Draft and chromosome-level genome assembly

A single D. avara individual from Girona, Spain (42.1145863 N, 3.168486 E) was used to generate a draft genome assembly. After extracting high molecular weight DNA, three PacBio Sequel II SMRT Cells (PacBio, Menlo Park, CA) were sequenced, generating a total of 32.1 Gb of data output in the form of CCS/HiFi reads (representing a coverage of ~55.8x based on the final assembly size). CCS/HiFi reads were assembled using the IPA pipeline, generating a draft assembly with a size of 648 Mb and a scaffold N50 of 3 Mb. The BUSCO score showed a completeness of metazoan single-copy orthologs of 86.8% and a duplication rate of 2.6%. By mapping the Omni-c libraries, generated by Cantana company, to the draft genome assembly, a high-quality chromosome-level genome assembly was generated for D. avara with a size 575 Mb, 162 scaffolds, and a scaffold N50 of 41Mb (Supplementary Table 3) structured in 15 chromosomes (Supplementary Figure 1). This assembly size is close to the k-mer based size prediction (499 Mb, Supplementary Figure 2). From the 147 scaffold sequences that were not assigned to a chromosome, 61 have a high similarity blast hit (E< 1e-30) to one of the assembled chromosomes, indicating they are alternative haplotype regions. From the remaining 86 sequences, after mapping against NCBI nr database using blast search, 37 sequences corresponded to microbial genomes, recovering the complete genomes of several bacteria, mostly of the genus Endozoicomonas (four different strains), and virus (Supplementary Table 4). The downstream analyses were carried out only with the sequences assigned to chromosome 1 to 15. The genome assembly was 64.9% complete according to the BUSCO score for metazoans (Figure 1B; Supplementary Table 3).

Figure 1

Figure 1. Dysidea avara species and its genome statistics. (A) Underwater photograph of the sponge D avara. (B) BUSCO scores indicating the quality of the genome and its annotation, taking into account the eukaryota_odb10 and metazoa_odb10 databases. The BUSCO analysis was run with Busco v.5 version.

The automated genomic annotation of this chromosome-level assembly identified 31,387 protein coding genes. We were able to assign KEGG pathways to 7,284 of these gene models, GO terms to 9,514 models and PFAM domains to 18,198 gene models (Supplementary Table 5). BUSCO results showed that our genome annotation strategy was able to retrieve up to 77.2% of single-copy metazoan orthologs (Figure 1B; Supplementary Table 3).

Comparing the Busco statistics of our genome with the ASG generated genome assembly of D. avara, we found that duplication rates of the genome assemblies were 2.3% and 2.0% for our genome and the ASG genome respectively (with completeness reaching 64.9% and 74.6% respectively) (Supplementary Table 3). When comparing the Busco results of the annotated genomes, the duplication level was 2.3% for our annotated genome while it was 44.3% for the ASG annotated genome (Supplementary Table 3).

3.2 Repertoire of protein domains related to animal immune receptors

We examined the presence of 24 conserved Pfam domains that are related to immune receptors in metazoans (21) and reported the number of genes that contain each domain (Table 1). We detected all the target domains in D. avara genome, except for the Leucine-rich repeat C-terminal domain, LRCNT (Table 1). The I-set (s) and V-set Pfam domains were represented by a relatively much higher number of gene sequences (637 sequences and 594 sequences, respectively) than other Pfam domains, such as the BIR, SEFIR, LRRNT and CTLD that were present in much less sequences (15 sequences, 1 sequence, 2 sequences, and 3 sequences, respectively). The C-set and NACHT domains were also represented by a relatively high number of sequences (371 and 230 sequences, respectively) (Table 1). Similar patterns were observed in most other sponge species of our study (Supplementary Table 6). All the sequences comprising LRR_2 and LRR_3 domains were already recovered in the LRR_1 domain search, while most of the LRR domains found in D. avara annotated genome corresponded to LRR_4 and LRR_6 domains, according to SMART search. We also found 14 gene sequences that contained TIR Pfam domain: 6 out of the 14 TIR domain-containing genes were combined with an immunoglobulin domain and one combined the TIR domain with a DEATH domain. Finally, 19 gene sequences contained TIR2 Pfam domain (Table 1), and 11 of those sequences were commonly assigned to both TIR and TIR2 domains.

Table 1

Table 1. Abundance of Pfam domains related to immunity in the genome of Dysidea avara.

3.3 Automated characterization of D. avara NOD-like receptors

The eggNOG automated annotation of D. avara genome identified 242 genes with sequence similarity to NLRs. Based on the domain architecture, 126 out of these 242 genes belonged to bona fide NLR genes (Supplementary Table 7). In particular, we identified 20 NLRC, 25 NLRD, and 81 NLRX genes (Figure 2A; Supplementary Table 7). NLRX includes mainly NACHT-LRR architectures without an identified N-terminal domain, but also a few cases of N-terminal domain different from CARD and DEATH (Figure 2A; Supplementary Table 7). Those 126 genes were spread homogeneously across the 15 chromosomes of D. avara genome, though some NLRX genes were found forming clusters in the chromosomes (e.g. in chromosome C8; Figure 2B; Supplementary Table 8). Finally, we identified another 63 genes that had only a NACHT annotation and another 64 NACHT containing genes with domain architectures not associated to bona fide NLRs (i.e. NACHT-WD40; NACHT-TPR; sushi-NACHT) (Supplementary Table 9); thus, we did not include those in our further analysis.

Figure 2

Figure 2. NOD-Like Receptors (NLRs) in Dysidea avara genome. (A) Domain architecture of proteins annotated to the different categories of bona fide NLRs (NLRD, NLRC, NLRX) in the genome assembly of D avara. The number indicates the annotated genes assigned to each NLR category. The graph was generated in Biorender.com. (B) Genes coding for bona fide NLRs were depicted in the chromosomes they are located in the genome. The graph was created with RIdeogram. NLRC, CARD-NACHT-LRR; NLRD, DEATH-NACHT-LRRs; and NLRX, NACH-LRR. CARD, Caspase recruitment domain; LRR, Leucine-rich repeat domain; Kix, kinase-inducible domain (KID) interacting domain; TNFR, tumor necrosis factor receptor; RPT, internal repeat; CASc, caspase; UBQ, ubiquitin; Litaf, LPS-induced TNF-activating factor; COR, C-terminal of Ras of Complex domain.

3.4 NLR categories in D. avara

We first focused on the NACHT containing proteins which make part of the different NLR categories in D. avara. We found that D. avara bona fide NLRs form two distinct large clusters separating the NLRC and NLRD categories (Figure 3). The NLRX category was present in both groups, but the majority of NLRX were more closely related to the NLRC category. Interestingly, those NLRX closer related to NLRC were clustered according to the number of C-terminal LRR repeats: a group of NLRX with low number of LRR repeats (2–4 repeats), another cluster with NLRX with 7–22 LRR repeats and finally a cluster with 15–30 LRR repeats (Figure 3).

Figure 3

Figure 3. Phylogenetic analysis of NLRs in Dysidea avara. The NACHT domain alignment from NLR sequences was used for the analysis (Supplementary File 1B; Supplementary File 2). Different clusters were shaped, depending on the NLR category NLRD or NLRC, while NLRX was grouped with the two other groups. Within the NLRX category, the number of C-terminal LRR repeats also determined different clusters. The ML tree was constructed with GTR Bootstrap expectation model and an estimated gamma shape parameter and 100 independent searches with RAxML. NLR, NOD-Like Receptor; LRR, Leucine Rich Repeat; NLRC, CARD-NACHT-LRR; NLRD, DEATH-NACHT-LRR; and NLRX, NACHT-LRR.

3.5 Phylogeny of poriferan bona fide NLRs

We then expanded the bona fide NLR search to the other 10 additional species, also including the recently available D. avara genome (Davar2, GCF_963678975.1) (Supplementary Table 1). Based on hmmer search with NACHT domain custom profile (Supplementary Table 2A; Supplementary File 1; Supplementary File 2A), we detected > 100 NACHT domain-containing proteins in each of the 11 sponges (Supplementary Figure 3B) but only a relatively small proportion (if any) of those proteins could be assigned to bona fide NLRs based on their domain architecture. For example, ca. 1/3 of NACHT domain-containing proteins in D. avara corresponded to bona fide NLRs, only one NLRX was detected in O. lobullaris, while no bona fide NLRs were detected at all in the hexactinellid sponge A. vastus, as so they were excluded from the downstream phylogenetic analysis (Supplementary Figure 3B). Dysidea avara had the largest genome span and the largest number of bona fide NLR genes (Supplementary Figure 3B). However, the number of bona fide NLR genes was not proportional to genome span in all species. In particular, the largest expansion relative to genome assembly span was found in H. panicea genome (Supplementary Figure 3). In all species with bona fide NLRs, most of them were assigned to the NLRX category, followed by NLRD category (Supplementary Figure 3). Looking into the phylogeny of the NACHT domain in bona fide NLR proteins across the different sponge species, we found that protein sequences clustered mostly by species and not by NLR category (Figure 4). Within the D. avara sequences (including those detected in both our genome and the genome generated by ASG), NLRC and NLRD were in different groups while the NLRX category appeared together with NLRC and NLRD, similar to what was described in Figure 3.

Figure 4

Figure 4. Phylogenetic analysis of NLRs from 8 sponge species based on the alignment of the NACHT domain. A phylogenetic analysis of NLRs was derived from the NACHT domain alignment in 8 sponge species, including our species (Supplementary File 1C; Supplementary File 2). As an outgroup, the NACHT domain from NLRs in the cnidarian Nematostella vectensis was used. The tree shows that NLRs are mainly grouped by species rather than by category. However, within the D. avara cluster, we see different clusters based on NLRD or NLRC category. The ML tree was constructed with GTR Bootstrap expectation model and an estimated gamma shape parameter and 100 independent searches with RAxML. NLRC, CARD–NACHT-LRR; NLRD, DEATH-NACHT-LRR; and NLRX, NACHT- LRR.

4 Discussion

This is the first study in which chromosome-level genomes of sponges were used to explore elements of innate immunity (D. avara genome only) and NLR expansions. Here, we generated the chromosome-level genome of the marine sponge D. avara and reported that this sponge contains almost the full diversity of conserved protein domains found in metazoan immune genes (Table 1), supporting the complex immune repertoire previously described in sponge draft genomes and transcriptomes (25, 26, 41, 64). Focusing further on NLR receptors, we observed an expansion of NACHT-containing genes in D. avara, detecting a large number of bona fide NLRs (NLRC, NLRD, NLRX) (Figure 3) as well as other NACHT domain-containing architectures previously described in ctenophora and cnidaria (23, 24, 65) (Supplementary Table 9). Finally, we reported a large expansion of NACHT-domain-containing genes in all sponges analyzed (Supplementary Figure 3). However, the bona fide NLR diversity was distinct for each species, with D. avara, A. oroides, H. panicea, and E. muelleri showing the largest bona fide NLR repertoires (Supplementary Figure 3; Figure 4).

4.1 Immune-related protein domains detected in D. avara

Dysidea avara chromosome-level genome comprises most of the protein domains characteristic of common metazoan immune receptors (Table 1), in agreement with reports in other sponge transcriptomes and genomes in this (Supplementary Table 6) and previous studies (26, 28, 40, 41). Some domains had larger gene expansions than others. I-set and V-set Pfam domains for instance were the most abundant in all sponge species studied here, but they were also one of the most abundant in ctenophores (65). The C-terminal LRRCT of the LRR structure was neither detected in our genome nor in any of the ASG annotated sponge genomes studied here (Table 1; Supplementary Table 6). LRRNT and LRRCT structures are mainly common in extracellular and membrane-associated LRR proteins (66). Interestingly, the ctenophores Mnemiopsis leidyi, Bolinopsis microptera, and Hormiphora californensis are also lacking the LRRCT sequence (65), indicating either that this part of the LRR structure appeared later in animal evolution or that there are different structure variants in early diverging metazoans.

The TIR domain is part of the Toll-like Receptors (TLRs) IL-1R, and Myd88 signaling gene (67–69). The canonical structure of the TLR receptor consists of the TIR domain combined with LRR as extracellular domains (16, 70, 71). However, this canonical structure is not present in sponges (26, 40, 72): the sponge homologous of the TIR domain in vertebrate TLRs, appears in combination with an extracellular immunoglobulin domain instead (40). According to that structure and homology, we identified 6 TLR-like genes in the D. avara genome, and the only TIR-DEATH sequence detected here is orthologous to vertebrate Myd88, according to the phylogenetic analysis conducted in a previous study (65). On the other hand, the TIR2 domain-containing genes of D. avara (Table 1), also found in the rest of studied sponge genomes (Supplementary Table 6) here and in previous studies (40) and previously also detected in Cnidaria (73), were homologous with those of Ctenophora, which have no signs of homology with TIR domains of TLR-like, or Myd88 structures of other metazoans (65). However, given that TIR2 domain-containing genes were overexpressed after challenging the ctenophore Mnemiopsis leidyi with pathogens (65), TIR2 containing proteins might also have similar functions in sponges as well.

4.2 An expanded Repertoire of NLRs is found in D. avara

The D. avara genome has a large expansion of NACHT domain-containing genes (126 bona fide NLR related genes and 128 additional genes) (Figure 3A; Supplementary Tables 7, 9). Such large expansions of this domain have previously been found in sponge draft genomes (25) and transcriptomes (28, 41, 64), ctenophores (65) and cnidarians (23, 24). Some of NACHT-containing genes have similar domain architectures in all early diverging metazoans (e.g. NACHT-WD40; NACHT-Zu5; NACHT-Ank) (Supplementary Table 9) (23, 24, 65) but still with unknown function. Though animal NLRs are reported as cytosolic (31), we confirmed here the presence of transmembrane domains in NLR architectures (Figure 2A), suggesting putative membrane-bound NLRs, as first reported in A. queenslandica (25). Bona fide NLR genes were distributed across D. avara chromosomes (Figure 2B), but we also observed genes of the same NLR category close together in the same chromosome (i.e. in chromosome 8, we saw several copies of NLRX nearby) (Figure 2B). Indeed, in 2 cases, we recovered almost identical genes close together in the same chromosomes (i.e. c6_44055680.417 with c6_44055680.432; 90% similarity; c13_30875229.36 with c13_30875229.44; 96% similarity; c2_56536150.10 with c2_56536150.14; 96.2% similarity), all belonging to NLRX category. Hence, we can hypothesize that NLR diversity in D. avara responds to gene duplication events. On another note, NLRX was also in close proximity with either NLRD or NLRC in most D. avara chromosomes (Figure 2B). This was also observed in the phylogenetic analysis of D. avara NLRs, in which NLRXs formed distinct clusters within either the NLRD or NLRC group (Figure 3). From the above, it is suggested that NLRD and NLC categories are true paralogues derived from NLRX. NLRXs have already been detected in Bacteria (as NB-ARC analogs), Fungi and protists (74) and they have been suggested as the ancestral core of NLR gene family, while the N-terminal DEATH/CARD was a later addition in the evolution of metazoan NLRs, providing the advantage of modulating signaling pathways of immunity and apoptosis (75).

4.3 NLR evolution in porifera

Our analysis of NLR diversity in additional chromosome-level sponge genomes revealed a high number of NACHT-domain containing genes in all species analyzed. The overall large expansion of immune gene families, including NLR-encoding genes, is a common characteristic of non-vertebrates (i.e. Hydra, sea anemone, sea urchin, amphioxus (76–79)), and this can be explained by the fact that they lack an adaptive immune system, and as such they rely on expansions of innate immune sensors to maintain effective immune defenses. In all cases, the NLRX category was the most expanded, followed by the NLRD category. Based on our pipeline, most of the species analyzed lacked NLRCs, except for D. avara (20 genes), H. panicea (1 gene), and C. crambe (1 gene). In a previous NLR analysis on the draft genome of A. queenslandica which included different gene models, Yuen et al. (25) identified 15 NLRC within a clade that contained NLRX and NLRD genes, with the latter group dominating A. queenslandica NLR diversity. In fact, the DEATH domain is part of apoptotic and immune regulators in basal metazoans, including sponges, and evolved earlier than the CARD domain, with the latter to be derived from the DEATH family (80). On the other hand, the CARD domain is more commonly found on the NLRs of bilaterians which make part of their inflammasome (81). This evolutionary trajectory agrees with the findings of a larger expansion of NLRD vs NLRC categories in our study (Figure 4; Supplementary Figure 3B; Supplementary Table 7), but also with reports in other marine invertebrates. Indeed, very few CARD domains were identified in the sea urchin genomes, and all the CARD-derived NLR functions of vertebrates were replaced by DEATH domain gene families, which dominate the genomes of sea urchin and amphioxus (76, 82). Under this evolutionary context, NLRs in Porifera, including NLRD and NLRC, were found in the most primitive position in the evolution of metazoan NLRs, which further expanded to either exclusively NLRD in certain invertebrates, or NLRC across most invertebrates and vertebrates (83). Overall, the presence of NLRC and NLRD already in sponge genomes suggests that these ancestral genes were already present in the last common ancestor of metazoans, making Porifera central to the origin of NLR evolution.

4.4 NLR repertoire in sponges is species-specific

The analysis of NLR categories revealed differences in the diversity of bona fide NLR repertoire among Porifera. The higher number of NACHT-domain containing genes and bona.

fide NLR genes recovered from the ASG D. avara genome in comparison to our genome.

(Supplementary Table 1; Supplementary Figure 3) are likely explained by the pipelines used.

for genome assembly and gene model predictions. Yet, both genomes support that D. avara has.

the largest number of bona fide NLRs (Supplementary Figure 3) and immune genes (Table 1; Supplementary Table 6) of all demosponges, pointing to a true biological difference. Following D. avara, A. oroides and H. panicea also presented large bona fide NLR expansions. Previously, Posadas et al. (64), studied the NLR repertoire in 16 sponge transcriptomes, identifying the largest number of NLRs in D. avara transcriptome (50 genes) as well. The variation in gene expansion and molecular divergence of NLRs among the different sponge species (Figure 4; Supplementary Figure 3) suggests different evolutionary trajectories to distinct selective pressures (84). The reduced repertoire of bona fide NLR compared to NACHT-domain containing genes may be explained by differences in the evolutionary history across species, divergent architectures of NACHT domain-containing genes and/or highly divergent LRR domains that could be missed in our analytical approach. Pathogen presence is a crucial selective agent in immune evolution, though here we analyzed species that occur in sympatry and at least in recent evolutionary times had likely encountered similar bacteria in the water while filter feeding. The maintenance of symbiotic host-microbe interactions can also dictate immune repertoires. Precisely, low microbial abundant (LMA) sponges, like D. avara, were suggested to have a larger expansion and higher diversity of NLR repertoire compared to HMAs (25, 85–87). For instance, the LMA sponge Stylissa carteri has a more expanded repertoire of immune domains, including the NACHT domain, than the HMA Xestospongia testudinaria (88). Large expansions in immune genes can also be correlated with variations in functions, or with allowing high expression levels (89, 90). Interestingly, some NLR genes were overexpressed while others were downregulated in the LMA D. avara while in contact with microbial-associated molecular patterns (a cocktail of bacterial lipopolysaccharide and peptidoglycan) (26) and non-symbiotic bacteria (42). On the other hand, NLRs did not participate at all during similar bacterial exposures in the HMA sponge Aplysina aerophoba (42). This means that NLRs could have different functions in different species and that the function may be related to the LMA-HMA status. However, we did not identify a clear pattern of HMA-LMA status on NLR expansion in our study (Supplementary Figure 3). Based on the species-specific trajectories suggested by our phylogenetic analysis, we propose that the evolution of NLR reflects species-specific traits, although functional convergence may rise from the LMA-HMA status.

5 Conclusion

The generation of a high-quality, chromosome-level genome for D. avara has enabled us to make a comprehensive investigation into the immune repertoire of this sponge, with a particular focus on NOD-like Receptors (NLRs). Our findings revealed a significant expansion of immune-related protein domains, including a notable presence of distinct NLR categories within D. avara. NLRX had the largest expansion, followed by NLRD, while the NLRC category was represented only by a few genes. Most likely, NLRs in sponges were diversified by NLRX. The species-specific grouping of NLRs across sponge genomes highlights that NLRs were diversified in sponges in order to respond to lineage-specific evolutionary patterns related to their immunity. This study provides valuable genomic resources, protocols, and novel insights into the early evolution of animal immune systems, reinforcing the importance of Porifera as a key group for understanding the origins and diversification of innate immunity.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

Ethics statement

The manuscript presents research on animals that do not require ethical approval for their study.

Author contributions

VK: Writing – original draft, Formal Analysis, Methodology, Writing – review & editing. MT-O: Writing – review & editing, Writing – original draft, Methodology, Formal Analysis. AM-G: Methodology, Writing – review & editing. AF: Writing – review & editing. JF: Methodology, Writing – review & editing. RS: Funding acquisition, Conceptualization, Writing – review & editing. UH: Supervision, Writing – review & editing. TR: Conceptualization, Funding acquisition, Supervision, Writing – review & editing. LP: Conceptualization, Supervision, Writing – review & editing, Funding acquisition, Formal Analysis, Writing – original draft.

Funding

The author(s) declared that financial support was received for this work and/or its publication. This work was funded by the Deutsche Forschungsgemeinschaft (DFG) Sequencing Call proposal IMMUBASE [417981041 to LP, RS and TR]; and by the DFG Collaborative Research Centre (CRC1182) “Origin and Function of Metaorganisms” (B2) [261376515 to RS and TR]. LP was supported by the “Ramón y Cajal” program (grant RYC2022-036761-I), funded by MICIU/AEI/10.13039/50/501100011033 and by ESF +. The authors acknowledge support of the publication fee by the CSIC Open Access Publication Support Initiative (PROA) through its Unit of Information Resources for Research (URICI).

Acknowledgments

We are grateful to. Dr. Ana Riesgo for her helpful support on tree reconstruction. We acknowledge the staff from CCGA sequencing facility in Kiel for the cDNA library preparation and sequencing as also their support.

Conflict of interest

The author(s) declared that this work was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2026.1725140/full#supplementary-material

References

1. Feuda R, Dohrmann M, Pett W, Philippe H, Rota-Stabelli O, Lartillot N, et al. Improved modeling of compositional heterogeneity supports sponges as sister to all other animals. Curr Biol. (2017) 27:3864–70. doi: 10.1016/j.cub.2017.11.008

PubMed Abstract | Crossref Full Text | Google Scholar

2. Reitner J and Wörheide G. Non-lithistid fossil demospongiae — Origins of their palaeobiodiversity and highlights in history of preservation. In: Hooper J, Van Soest R, and Willenz P, editors. Systema porifera: A guide to the classification of sponges. Springer US, Boston, MA (2002). p. 52–68.

Google Scholar

3. Love G, Grosjean E, Stalvies C, Fike D, Grotzinger J, Bradley A, et al. Fossil steroids record the appearance of Demospongiae during the Cryogenian period. Nature. (2009) 457:718–21. doi: 10.1038/nature07673

PubMed Abstract | Crossref Full Text | Google Scholar

4. Antcliffe J. The oldest compelling evidence for sponges is still early Cambrian in age – reply to Love and Summons (2015). Palaeontology. (2015) 58:1–3. doi: 10.1111/pala.12195

Crossref Full Text | Google Scholar

5. Van Soest RWM, Boury-Esnault N, Hooper JNA, Rützler K, de Voogd NJ, Alvarez B, et al. World porifera database. Available online at: http://www.marinespecies.org/porifera (Accessed October 16, 2020).

Google Scholar

6. Schmidt O. Grundzüge einer Spongien-Fauna des atlantischen Gebietes. Engelmann W, editor. Leipzig: Wilhelm Engelmann (1870), p. 42

Google Scholar

7. Sollas WJ. The sponge-fauna of Norway; a report on the rev. A.M. Norman’s collection of sponges from the norwegian coast. In: The Annals and magazine of natural history; zoology, botany, and geology, vol. 5. London: Taylor and Francis, Ltd (1880). p. 130–44.

Google Scholar

8. Bowerbank JS. On the anatomy and physiology of the spongiadae. Part II. Philos Trans R Soc. (1862) 152:747–829. doi: 10.1098/rstl.1862.0035

Crossref Full Text | Google Scholar

9. Bergquist PR. Sponges. Hutchinson. London: University of California Press, Berkeley & Los Angeles (1978) p. 1–268.

Google Scholar

10. Leys SP and Hill A. The physiology and molecular biology of sponge tissues. Adv Mar Biol. (2012) 62:1–56. doi: 10.1016/B978-0-12-394283-8.00001-1

PubMed Abstract | Crossref Full Text | Google Scholar

11. Hentschel U, Piel J, Degnan SM, and Taylor MW. Genomic insights into the marine sponge microbiome. Nat Rev Microbiol. (2012) 10:641–54. doi: 10.1038/nrmicro2839

PubMed Abstract | Crossref Full Text | Google Scholar

12. Webster NS and Taylor MW. Marine sponges and their microbial symbionts: Love and other relationships. Environ Microbiol. (2012) 14:335–46. doi: 10.1111/j.1462-2920.2011.02460.x

PubMed Abstract | Crossref Full Text | Google Scholar

13. Pita L, Erwin PM, Turon X, and López-Legentil S. Till death do us part: Sta ble sponge-bacteria associations under thermal and food shortage stresses. PloS One. (2013) 8:e80307. doi: 10.1371/journal.pone.0080307

PubMed Abstract | Crossref Full Text | Google Scholar

14. Webster NS and Thomas T. The sponge hologenome. mBio. (2016) 7:1–14. 8:e80307. doi: 10.1371/journal.pone.0080307

PubMed Abstract | Crossref Full Text | Google Scholar

15. Leys SP, Yahel G, Reidenbach MA, Tunnicliffe V, Shavit U, and Reiswig HM. The sponge pump: the role of current induced flow in the design of the sponge body plan. PloS One. (2011) 6:e27787. doi: 10.1371/journal.pone.0027787

PubMed Abstract | Crossref Full Text | Google Scholar

16. Dierking K and Pita L. Receptors mediating host-microbiota communication in the metaorganism: the invertebrate perspective. . Front Immunol. (2020) 11:1251. doi: 10.3389/fimmu.2020.01251

PubMed Abstract | Crossref Full Text | Google Scholar

17. Nyholm SV and Graf J. Knowing your friends: invertebrate innate immunity fosters beneficial bacterial symbioses. Nat Rev Microbiol. (2012) 10:815–27. doi: 10.1038/nrmicro2894

PubMed Abstract | Crossref Full Text | Google Scholar

18. Chu H and Mazmanian SK. Innate immune recognition of the microbiota promotes host-microbial symbiosis. Nat Immunol. (2013) 14:668–75. doi: 10.1038/ni.2635

PubMed Abstract | Crossref Full Text | Google Scholar

19. Bosch TCG. Rethinking the role of immunity: lessons from Hydra. Trends Immunol. (2014) 35:495–502. doi: 10.1016/j.it.2014.07.008

PubMed Abstract | Crossref Full Text | Google Scholar

20. Yang B, Yuen B, Degnan BM, and Degnan SM. An ancient role for conserved innate immunity transcription factors in holobiont self-recognition. bioRxiv. (2025) 2025:05.27.656460.

Google Scholar

21. Buckley KM and Rast JP. Characterizing immune receptors from new genome sequences. Methods Mol Biol. (2011) 748:273–98. doi: 10.1007/978-1-61779-139-0_19

PubMed Abstract | Crossref Full Text | Google Scholar

22. Rosenstiel P, Philipp EER, Schreiber S, and Bosch TCG. Evolution and function of innate immune receptors – insights from marine invertebrates. J Innate Immun. (2009) 1:291–300. doi: 10.1159/000211193

PubMed Abstract | Crossref Full Text | Google Scholar

23. Lange C, Hemmrich G, Klostermeier UC, Lopez-Quintero JA, Miller DJ, Rahn T, et al. Defining the origins of the NOD-like receptor system at the base of animal evolution. Mol Biol Evol. (2011) 28:1687–702. doi: 10.1093/molbev/msq349

PubMed Abstract | Crossref Full Text | Google Scholar

24. Hamada M, Shoguchi E, Shinzato C, Kawashima T, Miller DJ, and Satoh N. The complex NOD-like receptor repertoire of the coral acropora digitifera includes novel domain combinations. Mol Biol Evol. (2013) 30:167–76. doi: 10.1093/molbev/mss213

PubMed Abstract | Crossref Full Text | Google Scholar

25. Yuen B, Bayes JM, and Degnan SM. The characterization of sponge NLRs provides insight into the origin and evolution of this innate immune gene family in animals. Mol Biol Evol. (2014) 31:106–20. doi: 10.1093/molbev/mst174

PubMed Abstract | Crossref Full Text | Google Scholar

26. Pita L, Hoeppner MP, Ribes M, and Hentschel U. Differential expression of immune receptors in two marine sponges upon exposure to microbial-associated molecular patterns. Sci Rep. (2018) 8:16081. doi: 10.1038/s41598-018-34330-w

PubMed Abstract | Crossref Full Text | Google Scholar

27. Buckley KM and Rast JP. Diversity of animal immune receptors and the origins of recognition complexity in the deuterostomes. Dev Comp Immunol. (2015) 49:179–89. doi: 10.1016/j.dci.2014.10.013

PubMed Abstract | Crossref Full Text | Google Scholar

28. Degnan SM. The surprisingly complex immune gene repertoire of a simple sponge, exemplified by the NLR genes: A capacity for specificity? Dev Comp Immunol. (2015) 48:269–74. doi: 10.1016/j.dci.2014.07.012

PubMed Abstract | Crossref Full Text | Google Scholar

29. Schulenburg H and Ewbank JJ. The genetics of pathogen avoidance in Caenorhabditis elegans. Mol Microbiol. (2007) 66:563–70. doi: 10.1111/j.1365-2958.2007.05946.x

PubMed Abstract | Crossref Full Text | Google Scholar

30. Ting JPY, Lovering RC, Alnemri ES, Bertin J, Boss JM, Davis BK, et al. The NLR gene family: A standard nomenclature. Immunity. (2008) 28:285–7. doi: 10.1016/j.immuni.2008.02.005

PubMed Abstract | Crossref Full Text | Google Scholar

31. Proell M, Riedl SJ, Fritz JH, Rojas AM, and Schwarzenbacher R. The nod-like receptor (NLR) family: A tale of similarities and differences. PloS One. (2008) 3:e2119. doi: 10.1371/journal.pone.0002119

PubMed Abstract | Crossref Full Text | Google Scholar

32. Laroui H, Yan Y, Narui Y, Ingersoll SA, Ayyadurai S, Charania MA, et al. l-ala-γ-d-glu-meso-diaminopimelic acid (DAP) interacts directly with leucine-rich region domain of nucleotide-binding oligomerization domain 1, increasing phosphorylation activity of receptor-interacting serine/threonine-protein kinase 2 and its interaction with nucleotide-binding oligomerization domain 1. J Biol Chem. (2011) 286:31003–13. doi: 10.1074/jbc.M111.257501

PubMed Abstract | Crossref Full Text | Google Scholar

33. Mo J, Boyle JP, Howard CB, Monie TP, Davis BK, and Duncan JA. Pathogen sensing by nucleotide-binding oligomerization domain-containing protein 2 (NOD2) is mediated by direct binding to muramyl dipeptide and ATP. J Biol Chem. (2012) 287:23057–67. doi: 10.1074/jbc.M112.344283

PubMed Abstract | Crossref Full Text | Google Scholar

34. Philpott DJ, Sorbara MT, Robertson SJ, Croitoru K, and Girardin SE. NOD proteins: regulators of inflammation in health and disease. Nat Rev Immunol. (2014) 14:9–23. doi: 10.1038/nri3565

PubMed Abstract | Crossref Full Text | Google Scholar

35. Cuvillier-Hot V, Boidin-Wichlacz C, Slomianny C, Salzet M, and Tasiemski A. Characterization and immune function of two intracellular sensors, HmTLR1 and HmNLR, in the injured CNS of an invertebrate. Dev Comp Immunol. (2011) 35:214–26. doi: 10.1016/j.dci.2010.09.011

PubMed Abstract | Crossref Full Text | Google Scholar

36. Petnicki-Ocwieja T, Hrncir T, Liu YJ, Biswas A, Hudcovic T, Tlaskalova-Hogenova H, et al. Nod2 is required for the regulation of commensal microbiota in the intestine. Proc Natl Acad Sci. (2009) 106:15813–8. doi: 10.1073/pnas.0907722106

PubMed Abstract | Crossref Full Text | Google Scholar

37. Wiens M, Korzhev M, Perovic-Ottstadt S, Luthringer B, Brandt D, Klein S, et al. Toll-like receptors are part of the innate immune defense system of sponges (Demospongiae: porifera). Mol Biol Evol. (2007) 24:792–804. doi: 10.1093/molbev/msl208

PubMed Abstract | Crossref Full Text | Google Scholar

38. Srivastava M, Simakov O, Chapman J, Fahey B, Gauthier MEA, Mitros T, et al. The Amphimedon queenslandica genome and the evolution of animal complexity. Nature. (2010) 466:720–6. doi: 10.1038/nature09201

PubMed Abstract | Crossref Full Text | Google Scholar

39. Krishnan A, Dnyansagar R, Almén MS, Williams MJ, Fredriksson R, Manoj N, et al. The GPCR repertoire in the demosponge Amphimedon queenslandica: insights into the GPCR system at the early divergence of animals. BMC Evol Biol. (2014) 14:270. doi: 10.1186/s12862-014-0270-4

PubMed Abstract | Crossref Full Text | Google Scholar

40. Riesgo A, Farrar N, Windsor PJ, Giribet G, and Leys SP. The analysis of eight transcriptomes from all poriferan classes reveals surprising genetic complexity in sponges. Mol Biol Evol. (2014) 31:1102–20. doi: 10.1093/molbev/msu057

PubMed Abstract | Crossref Full Text | Google Scholar

41. Schmittmann L, Franzenburg S, and Pita L. Individuality in the immune repertoire and induced response of the sponge halichondria panicea. Front Immunol. (2021), 12:689051. doi: 10.3389/fimmu.2021.689051

PubMed Abstract | Crossref Full Text | Google Scholar

42. Marulanda-Gomez AM, Ribes M, Franzenburg S, Hentschel U, and Pita L. Transcriptomic responses of Mediterranean sponges upon encounter with symbiont microbial consortia. BMC Genomics. (2024) 25:674. doi: 10.1186/s12864-024-10548-z

PubMed Abstract | Crossref Full Text | Google Scholar

43. Schmidt O. Die Spongien des adriatischen Meeres. Leipzig. (1862) 1862:1–7.

Google Scholar

44. Uriz MJ, Rosell D, and Martín D. The sponge population of the cabrera archipelago (Balearic islands): characteristics, distribution, and abundance of the most representative species. Mar Ecology. (1992) 13:101–17. doi: 10.1111/j.1439-0485.1992.tb00343.x

Crossref Full Text | Google Scholar

45. Gloeckner V, Wehrl M, Moitinho-Silva L, Gernert C, Schupp P, Pawlik JR, et al. The HMA-LMA dichotomy revisited: an electron microscopical survey of 56 sponge species. Biol Bull. (2014) 227:78–88. doi: 10.1086/BBLv227n1p78

PubMed Abstract | Crossref Full Text | Google Scholar

46. Tommonaro G, Iodice C, AbdEl-Hady FK, Guerriero G, and Pejin B. The mediterranean sponge dysidea avara as a 40 year inspiration of marine natural product chemists. J Homeopathy Ayurvedic Med. (2014) s1:001. doi: 10.4172/2167-1206.S1-001

Crossref Full Text | Google Scholar

47. De Caralt S, Bry D, Bontemps N, Turon X, Uriz MJ, and Banaigs B. Sources of secondary metabolite variation in dysidea avara (Porifera: demospongiae): the importance of having good neighbors. Mar Drugs. (2013) 11:489–503. doi: 10.3390/md11020489

PubMed Abstract | Crossref Full Text | Google Scholar

48. Ribes M, Jiménez E, Yahel G, López-Sendino P, Diez B, Massana R, et al. Functional convergence of microbes associated with temperate marine sponges. Environ Microbiol [Internet]. (2012) 14:1224–39. doi: 10.1111/j.1462-2920.2012.02701.x

PubMed Abstract | Crossref Full Text | Google Scholar

49. Rix L, Ribes M, Coma R, Jahn MT, de Goeij JM, van Oevelen D, et al. Heterotrophy in the earliest gut: a single-cell view of heterotrophic carbon and nitrogen assimilation in sponge-microbe symbioses. ISME J. (2020) 14, 2554–2567. doi: 10.1038/s41396-020-0706-3

PubMed Abstract | Crossref Full Text | Google Scholar

50. Morganti TM, Ribes M, Yahel G, and Coma R. Size is the major determinant of pumping rates in marine sponges. Front Physiol. (2019) 10:1474. doi: 10.3389/fphys.2019.01474

PubMed Abstract | Crossref Full Text | Google Scholar

51. de Caralt S, Otjens H, Uriz MJ, and Wijffels RH. Cultivation of sponge larvae: settlement, survival, and growth of juveniles. Mar Biotech. (2007) 9:592–605. doi: 10.1007/s10126-007-9013-5

PubMed Abstract | Crossref Full Text | Google Scholar

52. Pita L, Fraune S, and Hentschel U. Emerging sponge models of animal-microbe symbioses. Front Microbiol. (2016) 7. doi: 10.3389/fmicb.2016.02102

PubMed Abstract | Crossref Full Text | Google Scholar

53. Wehrl M, Steinert M, and Hentschel U. Bacterial uptake by the marine sponge aplysina aerophoba. Microb Ecol. (2007) 53:355–65. doi: 10.1007/s00248-006-9090-4

PubMed Abstract | Crossref Full Text | Google Scholar

54. Putnam NH, Connell BO, Stites JC, Rice BJ, Hartley PD, Sugnet CW, et al. Chromosome-scale shotgun assembly using an in vitro method for long-range linkage. Genome Res. (2016) 26:342–50. doi: 10.1101/gr.193474.115

PubMed Abstract | Crossref Full Text | Google Scholar

55. Simão FA, Waterhouse RM, Ioannidis P, Kriventseva EV, and Zdobnov EM. BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs. Bioinformatics. (2015) 31:3210–2. doi: 10.1093/bioinformatics/btv351

PubMed Abstract | Crossref Full Text | Google Scholar

56. Vurture GW, Sedlazeck FJ, Nattestad M, Underwood CJ, Fang H, Gurtowski J, et al. GenomeScope: fast reference-free genome profiling from short reads. Bioinformatics. (2017) 33:2202–4. doi: 10.1093/bioinformatics/btx153

PubMed Abstract | Crossref Full Text | Google Scholar

57. Di Tommaso P, Chatzou M, Floden EW, Barja PP, Palumbo E, and Notredame C. Nextflow enables reproducible computational workflows. Nat Biotechnol. (2017) 35:316–9. doi: 10.1038/nbt.3820

PubMed Abstract | Crossref Full Text | Google Scholar

58. Letunic I, Khedkar S, and Bork P. SMART: recent updates, new developments and status in 2020. Nucleic Acids Res. (2021) 49:D458–60. doi: 10.1093/nar/gkaa937

PubMed Abstract | Crossref Full Text | Google Scholar

59. Hao Z, Lv D, Ge Y, Shi J, Weijers D, Yu G, et al. RIdeogram : drawing SVG graphics to visualize and map genome-wide data on the idiograms. PeerJ Comput Sci. (2020) 6:e251. doi: 10.7717/peerj-cs.251

PubMed Abstract | Crossref Full Text | Google Scholar

60. Edgar RC. MUSCLE: multiple sequence alignment with high accuracy and high throughput. Nucleic Acids Res. (2004) 32:1792–7. doi: 10.1093/nar/gkh340

PubMed Abstract | Crossref Full Text | Google Scholar

61. Kearse M, Moir R, Wilson A, Stones-Havas S, Cheung M, Sturrock S, et al. Geneious Basic: An integrated and extendable desktop software platform for the organization and analysis of sequence data. Bioinformatics. (2012) 28:1647–9. doi: 10.1093/bioinformatics/bts199

PubMed Abstract | Crossref Full Text | Google Scholar

62. Stamatakis A. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics. (2014) 30:1312–3. doi: 10.1093/bioinformatics/btu033

PubMed Abstract | Crossref Full Text | Google Scholar

63. Edler D, Klein J, Antonelli A, and Silvestro D. raxmlGUI 2.0: A graphical interface and toolkit for phylogenetic analyses using RAxML. Methods Ecol Evol. (2021) 12:373–7. doi: 10.1111/2041-210X.13512

Crossref Full Text | Google Scholar

64. Posadas N, Baquiran JIP, Nada MAL, Kelly M, and Conaco C. Microbiome diversity and host immune functions influence survivorship of sponge holobionts under future ocean conditions. ISME J. (2022) 16:58–67. doi: 10.1038/s41396-021-01050-5

PubMed Abstract | Crossref Full Text | Google Scholar

65. Koutsouveli V, Torres-Oliva M, Bayer T, Fuß J, Grossschmidt N, Marulanda-Gomez AM, et al. The chromosome-level genome of the Ctenophore Mnemiopsis leidyi A. Agassiz, 1865 reveals a unique immune gene repertoire. Genome Biol Evol. (2025) 17:evaf006. doi: 10.1093/gbe/evaf006

PubMed Abstract | Crossref Full Text | Google Scholar

66. Bella J, Hindle KL, McEwan PA, and Lovell SC. The leucine-rich repeat structure. Cell Mol Life Sci. (2008) 65:2307–33. doi: 10.1007/s00018-008-8019-0

PubMed Abstract | Crossref Full Text | Google Scholar

67. Bowie A and O’Neill LAJ. The interleukin-1 receptor/Toll-like receptor superfamily: signal generators for pro-inflammatory interleukins and microbial products. J Leukoc Biol. (2000) 67:508–14. doi: 10.1002/jlb.67.4.508

PubMed Abstract | Crossref Full Text | Google Scholar

68. Dunne A, Ejdebäck M, Ludidi PL, O’Neill LAJ, and Gay NJ. Structural complementarity of toll/interleukin-1 receptor domains in toll-like receptors and the adaptors mal and myD88. J Biol Chem. (2003) 278:41443–51. doi: 10.1074/jbc.M301742200

PubMed Abstract | Crossref Full Text | Google Scholar

69. Leulier F and Lemaitre B. Toll-like receptors — taking an evolutionary approach. Nat Rev Genet. (2008) 9:165–78. doi: 10.1038/nrg2303

PubMed Abstract | Crossref Full Text | Google Scholar

70. Takeda K, Kaisho T, and Akira S. Toll-like receptors. Annu Rev Immunol. (2003) 21:335–76. doi: 10.1146/annurev.immunol.21.120601.141126

PubMed Abstract | Crossref Full Text | Google Scholar

71. Coscia MR, Giacomelli S, Oreste U, and Oreste U. Toll-like receptors: an overview from invertebrates to vertebrates. ISJ. (2011) 8:210–26.

Google Scholar

72. Gauthier MEA, Du Pasquier L, and Degnan BM. The genome of the sponge Amphimedon queenslandica provides new perspectives into the origin of Toll-like and interleukin 1 receptor pathways. Evol Dev. (2010) 12:519–33. doi: 10.1111/j.1525-142X.2010.00436.x

PubMed Abstract | Crossref Full Text | Google Scholar

73. Poole AZ and Weis VM. TIR-domain-containing protein repertoire of nine anthozoan species reveals coral–specific expansions and uncharacterized proteins. Dev Comp Immunol. (2014) 46:480–8. doi: 10.1016/j.dci.2014.06.002

PubMed Abstract | Crossref Full Text | Google Scholar

74. Dyrka W, Lamacchia M, Durrens P, Kobe B, Daskalov A, Paoletti M, et al. Diversity and variability of NOD-like receptors in fungi. Genome Biol Evol. (2014) 6:3137–58. doi: 10.1093/gbe/evu251

PubMed Abstract | Crossref Full Text | Google Scholar

75. Koonin EV and Aravind L. Origin and evolution of eukaryotic apoptosis: the bacterial connection. Cell Death Differ. (2002) 9:394–404. doi: 10.1038/sj.cdd.4400991

PubMed Abstract | Crossref Full Text | Google Scholar

76. Huang S, Yuan S, Guo L, Yu Y, Li J, Wu T, et al. Genomic analysis of the immune gene repertoire of amphioxus reveals extraordinary innate complexity and diversity. Genome Res. (2008) 18:1112–26. doi: 10.1101/gr.069674.107

PubMed Abstract | Crossref Full Text | Google Scholar

77. Sodergren E, Weinstock GM, Davidson EH, Cameron RA, Gibbs RA, Angerer RC, et al. The genome of the sea urchin strongylocentrotus purpuratus. Sci (1979). (2006) 314:941–52. doi: 10.1126/science.1133609

PubMed Abstract | Crossref Full Text | Google Scholar

78. Putnam NH, Srivastava M, Hellsten U, Dirks B, Chapman J, Salamov A, et al. Sea anemone genome reveals ancestral eumetazoan gene repertoire and genomic organization. Sci (1979). (2007) 317:86–94.

PubMed Abstract | Google Scholar

79. Chapman JA, Kirkness EF, Simakov O, Hampson SE, Mitros T, Weinmaier T, et al. The dynamic genome of Hydra. Nature. (2010) 464:592–6. doi: 10.1038/nature08830

PubMed Abstract | Crossref Full Text | Google Scholar

80. Kohl A and Grütter MG. Fire and death: the pyrin domain joins the death-domain superfamily. C R Biol. (2004) 327:1077–86. doi: 10.1016/j.crvi.2004.08.006

PubMed Abstract | Crossref Full Text | Google Scholar

81. Franchi L, Warner N, Viani K, and Nuñez G. Function of Nod-like receptors in microbial recognition and host defense. Immunol Rev. (2009) 227:106–28. doi: 10.1111/j.1600-065X.2008.00734.x

PubMed Abstract | Crossref Full Text | Google Scholar

82. Hibino T, Loza-Coll M, Messier C, Majeske AJ, Cohen AH, Terwilliger DP, et al. The immune gene repertoire encoded in the purple sea urchin genome. Dev Biol. (2006) 300:349–65. doi: 10.1016/j.ydbio.2006.08.065

PubMed Abstract | Crossref Full Text | Google Scholar

83. Sun J, Leng J, and Song L. The evolution of NLR inflammasome and its mediated pyroptosis in metazoa. Int J Mol Sci. (2024) 25:11167. doi: 10.3390/ijms252011167

PubMed Abstract | Crossref Full Text | Google Scholar

84. Messier-Solek C, Buckley KM, and Rast JP. Highly diversified innate receptor systems and new forms of animal immunity. Semin Immunol. (2010) 22:39–47. doi: 10.1016/j.smim.2009.11.007

PubMed Abstract | Crossref Full Text | Google Scholar

85. Ryu T, Seridi L, Moitinho-Silva L, Oates M, Liew YJ, Mavromatis C, et al. Hologenome analysis of two marine sponges with different microbiomes. BMC Genomics. (2016) 17:1–11. doi: 10.1186/s12864-016-2501-0

PubMed Abstract | Crossref Full Text | Google Scholar

86. Germer J, Cerveau N, and Jackson DJ. The holo-transcriptome of a calcified early branching metazoan. Front Mar Sci. (2017) 4. doi: 10.3389/fmars.2017.00081

Crossref Full Text | Google Scholar

87. Schmittmann L, Rahn T, Busch K, Fraune S, Pita L, and Hentschel U. Stability of a dominant sponge-symbiont in spite of antibiotic-induced microbiome disturbance. Environ Microbiol. (2022) 24:6392–410. doi: 10.1111/1462-2920.16249

PubMed Abstract | Crossref Full Text | Google Scholar

88. Ryu T, Seridi L, Moitinho-Silva L, Oates M, Liew YJ, Mavromatis C, et al. Hologenome analysis of two marine sponges with different microbiomes. BMC Genomics. (2016) 17:158. doi: 10.1186/s12864-016-2501-0

PubMed Abstract | Crossref Full Text | Google Scholar

89. Sackton TB, Lazzaro BP, and Clark AG. Rapid expansion of immune-related gene families in the house fly, Musca domestica. Mol Biol Evol. (2017) 34:857–872. doi: 10.1093/molbev/msw285

PubMed Abstract | Crossref Full Text | Google Scholar

90. Zhang L, Li L, Guo X, Litman GW, Dishaw LJ, and Zhang G. Massive expansion and functional divergence of innate immune genes in a protostome. Sci Rep. (2015) 5:8693. doi: 10.1038/srep08693

PubMed Abstract | Crossref Full Text | Google Scholar