Brachiopod genome unveils the evolution of BMP signalling in bilaterian body patterning

Domain: biology
Published: 2026-03-12
Editorial outcome: major_revision

Abstract
The molecular control and ancestral state of dorsal–ventral patterning within spiralians remain unclear due to the remarkable diversity across species. Here we present a chromosome-level genome for the brachiopod Lingula anatina and apply functional transcriptomics to study dorsal–ventral patterning under BMP signalling control. We uncover asymmetrical activation of BMP signalling at the dorsal side of the gastrula, governed by ventral chordin expression and a ‘seesaw’ of BMP ligands. Using small-molecule drugs and recombinant proteins, we show that high BMP activity inhibits genes typically associated with neural patterning during gastrula and larval stages, similar to deuterostomes and non-spiralian protostomes. Our findings suggest deep conservation of this mechanism across all three major bilaterian clades, supported by striking similarities in BMP-regulated gene sets between brachiopods and Xenopus. We argue that the spiralian ancestor retained the ancestral bilaterian mechanism, although downstream network components have undergone developmental system drift.

Introduction
The integration of bone morphogenetic protein (BMP) signalling into the patterning of the dorsal–ventral axis is a cornerstone of bilaterian embryonic development. The last common ancestor of bilaterians likely established this axis through an interplay between BMP signalling molecules and their antagonists, such as Chordin1–3. In chordates, this interaction creates a BMP gradient, dividing the ectoderm into two distinct domains, a dorsal neural domain (low BMP) and a ventral epidermal domain (high BMP), in a process known as neural induction4,5. In the classical vertebrate model Xenopus, this gradient is regulated by the Spemann–Mangold organizer, a signalling centre that produces a cocktail of BMP antagonists, including Chordin4,6,7 and Noggin8,9. In the fruit fly Drosophila and other invertebrates, the BMP gradient is completely inverted, with dorsal Dpp/Bmp and ventral Sog/Chordin, but neural domain formation is still induced at the low BMP end of the gradient10–13. This reversed pattern is considered to represent the bilaterian ancestral state, shared by protostomes and ambulacrarian deuterostomes14–16—echinoderms and hemichordates—with a 180 degree axis inversion occurring in chordates15,17.

The action of several BMP ligands, including Bmp2/4, Bmp5–8 and Admp (anti-dorsalizing morphogenetic protein), combines to regulate the process of dorsal–ventral patterning. These proteins function as homodimers but more potently as heterodimers18–22, binding to type I and type II serine/threonine kinase receptors and resulting in the phosphorylation of the receptor-regulated Smad1/523–26. Nuclear-localized, phosphorylated Smad1/5 then binds with the common mediator Smad4 to act as a transcription factor, executing a transcriptional response to the BMP signal27,28. Notably, it is the expression domain of the antagonist Chordin, rather than those of the BMP ligands themselves, that is considered to be the primary determinant of the signalling gradient16,29. By binding BMP ligands and preventing their interaction with cell surface receptors, Chordin effectively attenuates BMP signalling29–32. While chordin is expressed at the low end of the BMP gradient, BMP ligand distribution varies significantly by taxon16, with bmp2/4 typically expressed at the opposite end to chordin3. In Drosophila, where BMP signalling is maximal dorsally, dpp, the ortholog of bmp2/4, is restricted to the dorsal side33, while screw, the ortholog of bmp5-8, is expressed ubiquitously34. In contrast, several distantly related animals, including Xenopus, the sea urchin Paracentrotus and even the sea anemone Nematostella, exhibit what has been denoted as a BMP ‘seesaw’, with BMP ligands produced at both the high and low ends of the signalling gradient35. In Xenopus, where BMP signalling is maximal ventrally, Bmp4 and Bmp7 are expressed ventrally but Bmp2 and Admp are dorsal5,35–37. In Paracentrotus, where BMP signalling is maximal dorsally, bmp2/4 and admp1 are coexpressed together with chordin on the ventral side, while admp2 is expressed dorsally16,38. In Nematostella’s directive axis, which may or may not be homologous to the dorsal–ventral axis of bilaterians29,39–43, the BMP ligands dpp and bmp5–8 are also co-expressed with chordin at the low BMP signalling end, with another BMP ligand, gdf5-like, produced at the opposite end39–41. The bimodal production of BMP ligands and associated negative feedback mechanisms, known as the BMP seesaw, is proposed to facilitate the production of a robust and precise gradient via self-regulatory mechanisms35. However, whether this arrangement was present in the ancestor of bilaterians or has evolved several times independently is currently uncertain29.

Aside from the dorsal–ventral axis inversion in chordates and the variability in BMP ligand expression, the ancestral BMP-mediated patterning mechanism remains highly conserved across diverse members of two of the three main groups of bilaterians: deuterostomes and ecdysozoans44,45. However, within Spiralia, a diverse group of protostomes including brachiopods, molluscs, annelids and platyhelminths that constitutes the third major bilaterian clade, there appears to be no consistent role for BMP signalling in dorsal–ventral axis patterning46. Spiralians display diverse patterning systems, including region-specific BMP signalling47–49 and even complete absence of dorsal–ventral BMP signalling50. For instance, perturbations to BMP signalling in the annelid Capitella teleta did not affect dorsal–ventral axis formation or central nervous system development but resulted in abnormal left–right asymmetries48,51. In the mollusc Crepidula fornicata, BMP signal perturbations affected neural and epidermal cell fates in the head region but not the trunk47. In another annelid, Chaetopterus pergamentaceus, Activin/Nodal signalling but not BMP signalling is the key specifier of the dorsal–ventral axis50,52. Indeed, the key BMP antagonist Chordin has been lost many times independently within annelids53 and other spiralians like platyhelminths54. In some species, such as the molluscs Lottia peitaihoensis (also known as Lottia goshimai) and Ilyanassa obsoleta, BMP signals involved in dorsal–ventral patterning are present but have a positive effect on neurogenesis55,56, opposite to that in ecdysozoans and deuterostomes. In the early-branching annelid Owenia fusiformis, BMP signals do not affect neurogenesis57, suggesting potential independent effects of BMP signals on dorsal–ventral fate and neural induction. In light of these experiments, the ancestral function of the BMP signalling pathway in spiralian development remains unclear51.

In this context, the focus of most spiralian studies on only two phyla, molluscs and annelids, limits our understanding. Broader sampling across spiralian phyla is essential to clarify the ancestral role of BMP signalling within this group and beyond. To this end, we aimed to study dorsal–ventral patterning in brachiopods, a phylum of marine spiralians with calcified dorsal–ventral shells that is relatively understudied compared to other spiralian phyla. Though previous studies have indicated that the brachiopods Novocrania anomala and Terebratalia transversa may use a BMP–Chordin network for dorsal–ventral patterning in a similar manner to ecdysozoans and deuterostomes58, understanding of this process remains incomplete. As a result, this group of animals makes ideal subjects for comprehensive studies of the BMP–Chordin axis to gain broader understanding of evolutionary scenarios regarding the ancestral role of BMP signalling and dorsal–ventral patterning within spiralians.

Here we present a chromosome-level genome for the brachiopod Lingula anatina. We then investigate dorsal–ventral patterning during L. anatina embryonic development using a combination of small-molecule drugs, recombinant proteins, functional transcriptomics and in situ hybridization. Specifically, we (i) characterize the BMP pathway gene repertoire; (ii) map the temporal and spatial expression of BMP ligands, antagonists and modulators during embryogenesis; (iii) identify the sites of BMP signal activation by these expression gradients; (iv) elucidate the role of asymmetric BMP signal activation in dorsal–ventral axis patterning; and (v) determine how BMP signalling influences neural domain formation.

Results
Chromosome-level assembly of the brachiopod genome

To explore the functions of BMP signalling in the brachiopod L. anatina, we first sequenced its genome to chromosome level using PacBio HiFi (405-fold coverage) and Hi-C (171 million paired reads) technologies (Supplementary Fig. 1, Supplementary Data 1, 2). The resulting assembly is 329.3 Mb in length with a 36.5% GC content and consists of 31.7% repetitive elements (Supplementary Fig. 2, Supplementary Data 3, 4). There have been recent expansions of repeats, particularly long terminal repeat elements, in the genome of L. anatina and other Lophophorata members (Supplementary Fig. 2). The assembly has a scaffold N50 of 30.1 Mb, a substantial improvement over the previous scaffold-level genome with a scaffold N50 of 0.5 Mb59, and is comparable in quality to other high-quality spiralian genomes53,60–62 (Fig. 1a, Supplementary Data 5).

The majority of the assembly consists of 10 chromosome-scale scaffolds accounting for 97.8% of the total length, consistent with previous karyotype observations in L. anatina (n = 10)63. Gene prediction and annotation using newly generated RNA sequencing data resulted in highly complete gene models with a BUSCO score of 98.1% (Metazoa odb10) (Supplementary Fig. 1, Supplementary Data 6–10). These gene models comprise 29,458 transcripts corresponding to 24,330 unique protein-coding genes, with extensive redundancy reduction achieved compared to the previous assembly (Supplementary Fig. 3, Supplementary Data 3). Gene family analysis reveals a gene content that is more stable than those of other Lophotrochozoans (Supplementary Fig. 4, Supplementary Data 11). Notably, we found that longer chromosomes have a higher density of protein-coding genes and a lower density of repeats (Supplementary Fig. 5). Phylogeny construction using exclusively genomic data (Supplementary Data 12) shows the brachiopod to be closely related to phoronids and bryozoans, as previously suggested64–71 (Fig. 1b). This result is robust to the inclusion of transcriptome data for the entoprocts Barentsia gracilis and Loxosomella nordgaardi72 and the cycliophoran Symbion pandora73 (Supplementary Fig. 6).

BMP pathway components in spiralian genomes

The absence of high-quality genomes has hindered accurate assessment of gene duplications and losses in numerous spiralian groups, especially in understudied lineages such as brachiopods, phoronids, bryozoans and nemerteans. To reconstruct the evolution of the BMP signalling pathway, we annotated BMP pathway components in L. anatina and thirteen additional spiralian genomes based on reported pathway genes in model systems such as Xenopus5, Drosophila74 and sea urchins75 and sought to characterize the evolution of BMP components through gene duplication and loss analysis (Fig. 1c, Supplementary Data 13, 14)76. In L. anatina, the BMP ligands bmp2/4, bmp3, bmp5-8 and admp were all present as single copies (Supplementary Fig. 7). Of the serine/threonine kinase receptors that receive BMP signals, bmpr1 and bmpr2 were single-copy while there was a partial duplication of acvr1 (Supplementary Fig. 8). Interestingly, the BMP signalling mediator smad1/5 was duplicated, while smad2/3, smad4 and smad6/7 were single-copy (Supplementary Fig. 9). We also annotated 11 BMP signalling modulators, including chordin, cv2 and bambi; of these, all were single-copy in L. anatina except noggin-like, which was not found (Supplementary Figs. 10–17). Thus, L. anatina shows a conserved BMP signalling toolkit, with most genes single-copy, just two duplications, and the loss of noggin-like, providing a reference point for comparison with other spiralians.

Characterising the BMP pathway gene sets across Spiralia, we found that the main ligands and modulators associated with dorsal–ventral patterning are highly conserved in the Lophotrochozoa, with bmp2/4, bmp5-8, admp and chordin present in all sampled species except specific annelid lineages53 (Fig. 1c). In contrast, other BMP components show remarkable variation considering their developmental importance. This is epitomised by bryozoans, which possess numerous lineage-specific duplications and losses to many key components at a scale comparable to known divergent groups like platyhelminths77. Two copies each of admp, cv2 and noggin are present and we failed to find acvr1, noggin-like and any of the three DAN family genes. We also found 11 copies of smad1/5 in the bryozoan Membranipora membranacea. Dynamic duplications and losses are not limited to bryozoans: the phoronid Phoronis australis has two copies of smad4 and three of cv2; the nemertean Lineus longissimus has two copies of bambi; two mollusc assemblies lack bmp3 and two annelid assemblies lack sclerostin. Interestingly, we identified nbl1 in L. anatina and P. australis but no other spiralian assemblies. Numerous gene duplicates are identified in the nemertean Notospermus geniculatus, but based on BUSCO duplication scores (Supplementary Fig. 1), the use of a low-quality draft assembly, and the lack of duplicates in the nemertean L. longissimus, these are likely attributable to false haplotypic duplications78,79. This finding highlights the need for high-quality chromosome-level genomes for gene content analyses. Overall, while L. anatina maintains a relatively conserved set of BMP pathway-related genes, repertoires exhibit notable evolutionary variation across other spiralians, indicating that BMP gradients may have evolved to accommodate diverse developmental modes and environmental contexts across various spiralian lineages.

Conserved syntenic associations of BMP pathway genes

We next questioned whether the developmental importance of BMP genes as an ancestral mechanism for body patterning in bilaterians has affected the evolution of their genomic position. To explore this hypothesis, we analysed genome-scale linkage conservation, focusing on the macrosyntenic positions of BMP genes in five bilaterians: the brachiopod L. anatina, the mollusc Pecten maximus61, the nemertean L. longissimus62, the annelid Owenia fusiformis53 and the chordate Branchiostoma floridae80. We first identified the correspondence between chromosomes and bilaterian ancestral linkage groups (ALGs)81, which are conserved to varying extents in spiralians82,83, then characterized the events of chromosome fusion and fission (Fig. 1d, Supplementary Data 15, 16). Strikingly, our analysis reveals extensive lineage-specific fusion-with-mixing in L. anatina, with chromosomes 1, 2, 3 and 7 all the products of fusion events (Fig. 1d). The giant 71.9 Mb L. anatina chromosome 1 consists of nine bilaterian ALGs (C1, G, I, J2, K, L, O1, O2 and R) and corresponds to six separate P. maximus chromosomes (Fig. 1d, Supplementary Fig. 18). Further comparison with the annelid O. fusiformis and nemertean L. longissimus suggests independent fusion events in brachiopods and annelids (Fig. 1d, Supplementary Fig. 18).

We next identified the genomic locations of BMP pathway components to assess their associations with an ALG. Our results show remarkable conservation: in each of the five selected bilaterian genomes, all 23 annotated BMP components are invariably linked to the same ALG (Fig. 1e, Supplementary Data 17). This exceeds the usual ALG conservation rate found in randomly selected genes (chi-square test; P = 0.003) (Supplementary Data 18), suggesting a higher tendency for BMP pathway genes to maintain their specific ALG associations. Even in species whose genomes have undergone extensive interchromosomal rearrangements, like that of the bryozoan M. membranacea84,85, 16 of 17 BMP pathway gene families present have conserved ALG associations. This suggests that BMP-related gene linkages are highly stable but not completely fixed, even in lineages where extensive interchromosomal rearrangements have occurred. This evolutionary conservation of the genomic locations of BMP pathway components may reflect the preservation of regulatory programs across large evolutionary spans.

To determine if the conservation of ALG associations is a characteristic unique to BMP pathway components, we extended our analysis to other genes with conserved, ancient roles in development. All 12 Wnt ligands found in both spiralians and chordates demonstrated 100% association with a single ALG across the five studied species, a rate significantly higher than the average ALG conservation (chi-square test; P = 0.038) (Supplementary Data 19). In addition, the Hox gene cluster is consistently associated with bilaterian ALG B2 in all analysed species. This persistent conservation of developmental genes' genomic locations, from brachiopods to chordates, despite widespread chromosomal rearrangements, implies evolutionary constraints on the mobility of these genes among ALGs. We propose that there is a stronger selective pressure at crucial developmental loci to maintain associations with regulatory elements in nearby genomic regions, thereby preventing their translocation and resulting in highly conserved macrosynteny.

Brachiopod embryos exhibit a BMP seesaw with dorsal BMP signal activation

The use of the BMP–Chordin signalling pathway for dorsal–ventral axis patterning is a conserved feature of the bilaterian body plan (Fig. 2a). This pathway operates through Smad1/5 phosphorylation, leading to nuclear translocation and regulation of downstream BMP signalling genes (Fig. 2b). To explore BMP–Chordin pathway gene expression during L. anatina development, we quantified the expression of BMP pathway genes in an RNA sequencing time course59 (Fig. 2c, d, Supplementary Fig. 19, Supplementary Data 20). Smad1/5 and bmp5–8 are expressed maternally, with their expression persisting through early development and into the larval stages (Fig. 2d). Although the L. anatina genome contains two copies of the smad1/5 gene (Supplementary Fig. 9), only one of these copies is active during early development (Fig. 2d), while the other is most highly expressed in the gonad (Supplementary Fig. 20). Core BMP pathway genes bmp2/4, bmp3 and chordin show early upregulation in the blastula, while modulators bambi and cv2 increase during the gastrula stage. The expression of admp peaks later, at the larval stage. BMP receptors are ubiquitously expressed, indicating that their expression pattern plays a non-crucial role in the timing of BMP pathway activation (Fig. 2). Interestingly, modulators like follistatin, gremlin, noggin and nbl1 show negligible expression in these stages, suggesting that they are not involved in dorsal–ventral patterning.

To elucidate the role of BMP signalling during L. anatina embryonic development, we first performed in situ hybridization for the four BMP signalling ligands (bmp2/4, bmp3, bmp5-8 and admp) and their main antagonist, chordin (Fig. 1e). Gene expression was assayed at three developmental stages: blastula, late gastrula and larva. In the blastula, chordin is expressed in the vegetal endomesoderm, bmp3 is localized to the dorsal ectoderm and bmp5-8 is expressed near-ubiquitously, but most highly in the dorsal ectoderm. Bmp2/4 is expressed to a low level in the endomesoderm. At the gastrula, there is a BMP seesaw of strong dorsal bmp2/4 and bmp3 expression and equally strong opposing ventral admp expression, in the same location as chordin. Expression of bmp2/4, bmp3 and bmp5-8 is also observed in the endomesoderm. In the larva, which is the product of significant reorganization post-gastrulation (Supplementary Fig. 21), bmp2/4 is expressed in the centre of the tentacle and the gut, while bmp3 and bmp5-8 are also expressed in the gut. Admp is expressed in the lateral part of the tentacle while chordin expression is near-absent but may persist at a low level in the same domain as admp.

We next investigated how these signalling ligand and antagonist expression domains translate to localized or asymmetrical activation of the BMP signalling pathway during L. anatina development. We assessed BMP pathway activity by detecting the presence or absence of phosphorylated Smad1/5 (pSmad1/5), which serves as a readout of BMP signalling. Using an antibody targeting a conserved PHNPISSVS phosphopeptide sequence present in metazoan Smad1/5 proteins (Supplementary Fig. 22), we examined the distribution of pSmad1/5 across a developmental time course. From the 1-cell to 128-cell stages, we detected no pSmad1/5 nuclear signal (Fig. 2d, f–n, Supplementary Fig. 23). Expression of the BMP ligand bmp2/4 is very low before the early blastula stage (Fig. 2d), suggesting that BMP signalling is not activated until after the cleavage stages. However, we did not directly test the effects of BMP signal manipulation at this stage, and further experimentation is required to confirm this. Nuclearized pSmad1/5 appears abruptly in the blastula, where it is asymmetrically localized at the animal pole towards the dorsal side (Fig. 2k). This occurrence aligns with the upregulation of bmp2/4, bmp3 and chordin (Fig. 2d) observed in RNA-seq data, suggesting an importance of these proteins to the regulation of this signal. The presence of pSmad1/5 then persists on the dorsal side through the gastrula and early larval stages (Fig. 2l–n), demonstrating dorsal activation of the BMP signalling pathway. This pattern is consistent with those found in other protostomes and ambulacrarian deuterostomes, indicating the presence of a deeply conserved dorsal–ventral asymmetry of BMP pathway activation in spiralians.

Asymmetric BMP activation regulates dorsal–ventral patterning

To investigate the role of asymmetric BMP activation, we perturbed BMP signalling with receptor inhibitors and recombinant BMP proteins to block or enhance its activity, respectively. Treatments were applied either from early blastula to late gastrula, when dorsal–ventral molecular features first appear, or extended to the one-pair-cirri larva, after axis patterning is established (Fig. 3a, Supplementary Data 21). In control late gastrulae, pSmad1/5 staining was restricted to the dorsal side (Fig. 3b,c). Inhibition of BMP receptors with LDN-193189 (LDN)86 or K02288 (K02)87, collectively referred to herein as BMP(-), abolished pSmad1/5 signal entirely (Fig. 3d–e). Conversely, treatment with recombinant mouse BMP4 protein (BMP)88, referred to as BMP(+), expanded nuclear pSmad1/5 across the whole gastrula (Fig. 3f), validating the efficacy of the manipulations. pSmad1/5 staining confirmed that both inhibitors and recombinant proteins remained effective until the larval stage (Supplementary Fig. 24). Notably, both BMP(-) and BMP(+) treatments disrupt gastrulation. BMP(-) treatment results in a shortened archenteron (Fig. 3d) or a complete delay in invagination when using a stronger BMP inhibitor (Fig. 3e). In contrast, BMP(+) treatment delays invagination but promotes ingression (Fig. 3f). These observations show the importance of the BMP gradient in gastrula development.

We next explored the influence of BMP signalling on cellular proliferation and the cell cycle by quantifying gastrula cells marked with the mitosis-specific phospho-histone H3. BMP(-) treatments significantly increased proliferative cell numbers compared to the control, whereas BMP(+) embryos showed a decrease in mitotic cells (P < 0.05) (Fig. 3g–j, o). This is consistent with dose-dependent effects of BMP signals reported from vertebrate development, where low doses promote proliferation while high doses inhibit proliferation and induce terminal differentiation89–91, and suggests that a similar effect may be exerted during L. anatina development.

During L. anatina development, the embryo develops a mantle fold which is subsequently split into two mantle lobes that go on to form the two shells (Fig. 3b)92. We assessed the necessity of BMP signalling to this process by manipulating the BMP pathway and observing embryonic shell formation with a chitin-binding probe. In control embryos, the larval shell is marked by chitin staining across both mantle lobes, excluding the dorsal edge (Fig. 3k). In contrast, BMP(-) embryos showed fused shell fields and a single, unseparated shell, indicated by continuous chitin staining, including at the dorsal edge (Fig. 3l–m). This suggests that the depressions that split the mantle fold into two lobes92 during normal development have failed to form properly. Meanwhile, BMP(+) embryos exhibited extended shell domains, with the chitinous shell extending ventrally, resulting in the tentacle and mouth being encased inside the shell (Fig. 3n). The lack of a clear opposite effect on shell formation in BMP(-) and BMP(+) embryos may be because the effects of BMP signalling on mantle lobe and shell formation are a product of interference with embryonic morphogenetic processes (e.g., the formation of the furrows that split the mantle fold) rather than a direct promotion or inhibition of shell formation. The disturbance of the distribution of proliferative cells may also contribute to this process. Overall, these experiments indicate that dorsal BMP signalling is essential for splitting the developing larval shell into two distinct valves59 and, more broadly, maintaining embryonic structural integrity and proper morphogenesis.

BMP signals control neural and cell cycle gene expression

We further explored the role of BMP signalling in L. anatina development through functional transcriptomics, using RNA sequencing on control, BMP(-) and BMP(+) embryos at both the late gastrula and one-pair-cirri larval stages (Supplementary Data 22–24). Differential gene expression and hierarchical clustering analyses identified gene sets regulated by BMP signalling. Genes downregulated under BMP(-) treatment and upregulated under BMP(+) treatment were classified as BMP-upregulated, whereas genes showing the opposite pattern were categorized as BMP-downregulated. At the late gastrula stage, this produced a set of 881 genes with strongly significant expression changes (fold-change > 4; P < 0.001) (Fig. 4a, Supplementary Fig. 25). Results were consistent across both independent BMP inhibitor molecules (LDN and K02), suggesting that the observed phenotypes were due to BMP signal disruption rather than off-target effects of a specific inhibitor.

We used gene ontology (GO) analysis to characterize these gene sets (Supplementary Data 25–28). Among BMP-downregulated genes, 'nervous system development' was one of the most highly enriched terms (Supplementary Fig. 25). To investigate this further, we sought to use a highly robust, conservative dataset of genes for which we had strong evidence for regulation by BMP signalling. To this end, we restricted the dataset to only genes that were BMP-regulated in both our late gastrula and early larval stage experiments, resulting in a total of 49 BMP-upregulated and 45 BMP-downregulated genes (Fig. 4b, c, Supplementary Fig. 26). In this dataset, the BMP-downregulated targets were dominated by genes demonstrating either exclusively neural expression or known to be expressed in neurons in specific developmental contexts in other model systems, including elav, which may be a universal neural marker93, foxb1, lhx1/5, netrin, nkx2.4, otx2, pax6, six3/6, sim1/2, soxb1 and soxb294–101 (Fig. 4h, Supplementary Figs. 27, 28).

To validate the transcriptomic data, we examined expression of four key neural markers—soxb1, soxb2, pax6 and elav—by in situ hybridization in control, BMP(+) and BMP(-) embryos (Fig. 4d, Supplementary Data 29). In control embryos, all four genes were expressed ventrally, directly opposite the dorsal BMP maximum (Fig. 4d, top row). For soxb1, soxb2 and pax6, inhibition of BMP signalling caused a dorsal expansion of expression domains, with stronger effects under K02 compared to LDN treatment (Fig. 4d, second and third rows). In contrast, BMP(+) treatment abolished detectable expression of these three genes (Fig. 4d, bottom row). Elav expression was slightly reduced under BMP(+), and showed a slight, though less clear, upregulation with K02 (Fig. 4c). This is consistent with the RNA-seq data, where the effect on elav is less pronounced than other neural genes. Overall, these data show that key neural genes are expressed in the same location as chordin, diametrically opposed to the high point of BMP signal activation. Moreover, their expression domains expand when BMP signalling is inhibited and are eliminated when it is overactivated, confirming that BMP signals restrict the expression of neural genes. This finding, supported by both transcriptomic and in situ hybridization data, contrasts with a recent study that reported an unusual positive effect of BMP signals on neural domain formation in spiralians55,56, but is consistent with the observed inhibitory effect of BMP in the neuroectoderm of non-spiralian protostomes102 and deuterostomes15,45. Our results suggest that inhibition of neural development by elevated BMP levels is a mechanism shared between spiralians, ecdysozoans and deuterostomes, despite the evolutionary diversification of embryological development modes in molluscs55 and annelids48.

In addition, at the gastrula stage, each of the 12 most statistically significant GO terms for BMP-downregulated genes relates to processes involved in cell cycle regulation (Supplementary Data 30), including ‘DNA replication’, ‘chromosome organization’ and ‘cell cycle’. BMP-downregulated genes include atr, cdc16, ctf18, pcna and pole, all of which encode essential components of the cellular DNA replication machinery103–108. All six mcm2–7 genes, encoding the core components of the eukaryotic replicative helicase complex109, are also BMP-downregulated. Conversely, no GO terms associated with BMP-upregulated genes relate to DNA replication or the cell cycle (Supplementary Data 25). This is consistent with the above result that the presence of mitotic cells was reduced in BMP(+) embryos and increased in BMP(-) embryos, and supports the contention that BMP signals inhibit cellular proliferation at this time point. A role for BMPs in regulating proliferation is well-known but the extent and direction of effect is highly context-dependent110.

Deep conservation of the dorsal–ventral patterning pathway across bilaterians

We next investigated whether the downstream effects of BMP signalling on gene expression that we observed in L. anatina show similarities to those of other bilaterians. We selected the vertebrate model Xenopus for this comparison, as the role of BMP signalling in dorsal–ventral patterning is well characterized5 and its phylogenetic position within the deuterostomes means that evolutionary inferences can be made that span the entire Bilateria. In Xenopus, a ventral signalling centre produces the BMP ligand Bmp4111 plus feedback inhibitors like Cv2112, Twisted-gastrulation113–115 and Bambi116, while the dorsal Spemann–Mangold organizer secretes Admp117 and a set of BMP antagonists that includes Chordin4,6,7, Follistatin118 and Noggin8,9. These antipodal signalling centres construct the dorsal–ventral BMP gradient that patterns the axis, with the neural inducers secreted from the Spemann–Mangold organizer causing the formation of neural tissue on the dorsal side.

We found that several brachiopod homologues of key Xenopus dorsal–ventral patterning genes were up- or down-regulated at the RNA level in the above BMP signal manipulation experiments. Specifically, BMP signalling in L. anatina downregulates homologues of the Xenopus Spemann–Mangold organizer genes chordin and admp, while upregulating homologues of the Xenopus ventral centre genes bambi and cv2 (Fig. 4b, c).

The parallel in regulatory logic between Xenopus and L. anatina suggests a conserved BMP signalling network architecture (Fig. 5a). We next investigated whether the spatial expression patterns of these genes in Xenopus are also conserved in L. anatina. In situ hybridization experiments show that bambi and cv2 are dorsally localized in L. anatina gastrulae, while chordin and admp are found ventrally (Fig. 5b). This is consistent with their expression in Xenopus when the axis inversion is taken into account. Notably, both chordin and admp were expressed at the boundary of the ventral ectoderm and endomesoderm, similar to the expression pattern observed in the cephalochordate B. floridae15 (Fig. 5b). Bambi and cv2 were expressed in domains of active BMP signals on the dorsal side. BMP(-) treatment causes expansion of chordin and admp expression domains and the near absence of bambi and cv2 expression. In contrast, BMP(+) treatment results in the near absence of expression of admp and chordin while bambi and cv2 are expressed universally (Fig. 5b). These experiments support our transcriptomic results showing that BMP signalling upregulates bambi and cv2 while downregulating chordin and admp. This result is also consistent with regulatory interactions observed in Xenopus, where BMP signalling represses chordin5 and admp117,119 while activating bambi116,120 and cv25,121. Overall, we found conservation of Xenopus dorsal–ventral patterning gene regulatory logic and gene expression patterns in the brachiopod L. anatina (Fig. 5c, d), indicating deep conservation of these pathways in spiralians.

Despite this remarkable conservation, several Xenopus Spemann–Mangold organizer and ventral centre genes (e.g., twsg, noggin and follistatin) do not appear in our restricted differentially expressed lists (Fig. 4b, c). We hypothesised that this may reflect either our highly conservative approach to differentially expressed gene identification or genuine divergence in the role of these genes between species. To investigate this, we studied the expression of 15 key Xenopus Spemann–Mangold organizer and ventral centre genes under each manipulation condition. Intriguingly, many of these genes are differentially expressed in response to BMP signalling in L. anatina, but not in the predicted direction. For instance, the Spemann–Mangold organizer gene gremlin is BMP-upregulated when it might be expected to be downregulated like chordin and admp, and vice versa for twsg (Fig. 5e). These findings, while requiring further validation, suggest that although the central BMP–Chordin axis is conserved, its downstream effects are evolutionarily variable. Accordingly, a phylostratigraphy approach using transcriptome age index (TAI) analysis122 shows that increasing the level of BMP signalling significantly decreases the evolutionary age of the transcriptome during both gastrula and larval stages of development (two-sided t test P < 0.05; Fig. 5f, g, Supplementary Data 31–33). This result suggests that activating the BMP pathway in L. anatina leads to the upregulation of evolutionarily young genes. Indeed, our findings show that during the developmental time course, the late blastula and early gastrula stages—following BMP signal activation—exhibit the highest TAI, indicating the youngest transcriptome (Supplementary Fig. 29). This aligns with transcriptome ages at similar stages in zebrafish development122. When combined with our BMP signal inhibition experiments, the findings suggest that BMP signals are involved in a gastrulation process that is marked by increased expression of genes that emerged relatively recently in evolution.

Methods
Biological materials and artificial fertilization

Between 2013 to 2016, gravid adults of the brachiopod L. anatina (approximately 150 to 200) were collected each August during low tide from Kasari Bay, Amami Island, Japan (28.440583 N, 129.667608 E). Around 50 individuals were housed in a 10 L aerated seawater tank and fed daily with 25 mL/tank of Chaetoceros calcitrans (5×107 cells/mL). Oocyte maturation was induced by injecting each gonad with 30 μL of 40 mM dibutyryl-cAMP dissolved in phosphate-buffered saline (PBS)59,63. To control fertilization timing, post-injection individuals were isolated in separate Petri dishes. Artificial spawning was facilitated by adjusting the temperature from 25°C to 29°C over 2 h and then rapidly returning it to 25°C for a cold shock144. Prior to inducing fertilization, sperm activity was inspected under a stereomicroscope. Fertilization efficacy decreased sharply 4 h after spawning. Depending on seasonal factors, spawning rates varied from 20% to 70%.

Genomic sample acquisition and DNA extraction

In August 2021, fresh adults of L. anatina were collected with the objective of obtaining high-quality genomic DNA suitable for both long-read and Hi-C sequencing methodologies. The dissected tissues, including the lophophore, mantle and adductor muscle, were immediately snap-frozen in liquid nitrogen. For genome resequencing, approximately 2 g of lophophore and mantle tissue and 1.5 g of adductor muscle tissue were collected from a single adult specimen. These tissues were then divided into two equal portions. One portion was used for genomic DNA extraction, while the other was reserved for Hi-C sequencing. The cetrimonium bromide (CTAB) method was utilised to extract high molecular weight genomic DNA. After extraction, the DNA was purified using the Blood & Cell Culture DNA Midi Kit (Qiagen 13343). DNA purity was assessed with a NanoDrop spectrophotometer, and DNA concentration was determined using the Qubit 4.0 Fluorometer. Genomic DNA quality was assessed by gel electrophoresis. Samples meeting the quality criteria showed a main band centred around 100 kb with no detectable signal below 40 kb. DNA purity was evaluated using NanoDrop (A260/280 = 1.8–2.0; A260/230 = 2.0–2.2), and DNA concentration was measured using both NanoDrop spectrophotometer and the Qubit 4.0 Fluorometer.

RNA sequencing of adult tissues

Total RNA was extracted from the mantle, gonad, lophophore, adductor muscle and pedicle using TRIzol. The RNA was reconstituted in 50 μL RNase-free water, and its concentration was determined with a Nanodrop spectrophotometer. RNA quality was verified using the 5400 Fragment Analyzer System (Agilent Technologies). Total RNA from each tissue sample was used directly for mRNA enrichment and library preparation using the TruSeq RNA Library Prep Kit v3, following the manufacturer’s protocol. Briefly, mRNA was purified using oligo(dT)-attached magnetic beads, followed by fragmentation and first-strand cDNA synthesis with random hexamer primers. Second-strand synthesis was then performed, and the resulting cDNA underwent end repair, A-tailing, adapter ligation, size selection, amplification and purification to complete library construction. Libraries were sequenced as 150 bp paired-end reads on the Illumina NovaSeq 6000. Quality control was conducted using FastQC145. The reads were trimmed with Trimmomatic (v0.39)146 and then aligned to the genome using STAR (v2.7.10b)147 for gene prediction. For gene expression analysis, transcript abundances were quantified using kallisto (v0.43.0)148 (Supplementary Data 34). The same method was used to process data from a published developmental dataset of L. anatina59 (Supplementary Data 35).

Genome sequencing and assembly

For PacBio HiFi circular consensus sequencing, SMRTbell libraries were constructed following PacBio's standard protocol, utilising the 15 kb preparation solutions. In brief, high-molecular-weight DNA was extracted from the lophophore sample and was used for DNA library preparations. The genomic DNA was sheared to the desired fragment size using g-TUBEs (Covaris). After removal of single-strand overhangs, damage repair, end-repair and A-tailing, DNA fragments were ligated with the hairpin adaptor. Post-ligation, the library was treated with nuclease, cleaned using the SMRTbell Enzyme Cleanup Kit, and purified with AMPure PB Beads (Beckman Coulter). The desired fragments were subsequently isolated using BluePippin (Sage Science). An Agilent 2100 Bioanalyzer was used to determine the size distribution of the library fragments. The final sequencing of the SMRTbell library was conducted on the PacBio Sequel II platform, utilising Sequencing Primer V2 and the Sequel II Binding 2.0 Kit at Grandomics Biosciences. HiFi reads with a mean length of 13 kb were assembled using Hifiasm (v0.16.0)149, and redundant haplotype sequences were removed with Purge_dups (v1.2.5)150.

The Hi-C library was constructed following a previously described method151. In brief, mantle, lophophore, and adductor muscle tissues were cut into 2 cm segments and cross-linked using a nuclei isolation buffer with 2% formaldehyde. Post cross-linking, tissues were ground to produce a nuclei suspension. The isolated nuclei were digested with 100 units of MboI and labelled with biotin-14-dATP. Unligated fragments had their biotin removed using T4 DNA polymerase. The DNA was sheared to 300–600 bp fragments, underwent blunt-end repair and A-tailing, and was isolated with streptavidin beads. After quality checks with a Qubit Fluorometer and an Agilent 2100 Bioanalyzer, the libraries were sequenced as 150 bp paired-end reads on the Illumina NovaSeq 6000 platform.

Hi-C data was used for assembly scaffolding152 using the following pipeline: trimmomatic (quality filtering)146; FastQC (quality control)145; bwa_mem2 (sequence alignment)153; Juicer (data filtering)154; 3D-DNA (Hi-C-assisted assembly)155; and JuiceBox (visualization and error correction)156. Having obtained a preliminary assembly, blasts, read coverage and GC content were used to check for sequences belonging to contaminants and symbionts, The assembly was first subjected to Diamond157 blastx against the nr database of NCBI with ‘sensitive mode’ parameter and an E-value of 1e-25. Subsequently, Hi-C reads were mapped to the original genome using HISAT2158 with default parameters and in single-end mode to calculate the coverage of each fragment. These three data types were visualised using blob plots created with BlobToolKit159. This methodology revealed a small cluster of short sequences with an unusually high GC content (greater than 45%) which blasted to non-brachiopod taxa. All such fragments were removed from the genome. A second blob plot, produced from the remaining fragments, identified only sequences blasting to Brachiopoda. Sequences present after this filtering process represent the final L. anatina assembly.

Gene prediction and annotation

Repetitive elements were annotated de novo using RepeatModeler (v2.0.4)160. RepeatModeler employs RepeatScout161 and RECON162 to identify transposable elements and also uses the long terminal repeat (LTR)-specific tools LTRharvest163 and LTR_retriever164. The BRAKER pipeline (v3.0.3)165–173 was then used for gene prediction and annotation. Genomes were first soft-masked with RepeatMasker (v4.1.5; sensitive mode)174, and then BRAKER was run using hints from mapped RNA sequencing data. RNA-seq reads (Supplementary Data 36) were downloaded from the NCBI Sequence Read Archive, trimmed with Trimmomatic (v0.39) and aligned with STAR (v2.7.10b) before input to BRAKER in BAM format. The RNA-seq datasets generated in this study were also used as hints for BRAKER. Gene annotation quality was assessed using BUSCO (v5.4.7)175. InterProScan176, KofamScan177 and EggNOG-mapper178,179 were used for functional annotation. For KofamScan, output is limited to hits where threshold > score (adjudged to be a significant hit). Orthologues of L. anatina genes in human and mollusc (Patella vulgata) genomes were identified with OrthoFinder (v2.5.4)180. RepeatLandscape was used to create Kimura substitution level plots for repeats in lophophorate genomes. Ribosomal RNA (rRNA) genes were annotated with barrnap (v0.9). Gene density and repeat density plots were made with RIdeogram (v0.2.2)181.

Phylogenomic analysis

Proteomes were downloaded from NCBI or produced by the gene prediction method outlined above (Supplementary Data 12). Orthology assignment was performed with OrthoFinder 175. OrthoSNAP (v0.01)182 was then run with default parameters to recover additional orthologues for phylogenetics. The resultant dataset contained 2,036 OrthoSNAP orthogroups. To determine the strength of the phylogenetic signal possessed by each orthogroup, we first aligned sequences with MAFFT (v7.520)183,184 and trimmed alignments with ClipKIT (v1.4.1)185. Individual gene trees were then constructed with IQ-TREE (v2.2.2.3)186 using ModelFinder automated model selection187 and UFBoot2 ultra-fast bootstrapping188. Orthogroups with an average bootstrap score over 85 % (n = 109) were selected for species tree building and alignments concatenated with PhyKIT (v1.11.7)189. The final tree was built using IQ-TREE as above (model LG+F+R5).

Comparative genomics

We used OrthoFinder to identify putative homologues of BMP pathway components in lophotrochozoan genomes. Genes were identified by orthology to sequences from human and Drosophila, which are well-annotated and verified. All cases of putative gene losses and duplications were then manually verified using gene tree construction with IQ-TREE (v2.2.2.3) plus reciprocal blast searches and microsynteny comparisons where necessary. Gene trees for each superfamily were built using a pipeline of MAFFT (v7.520) alignment, ClipKIT (v1.4.1) trimming and IQ-TREE (v2.2.2.3) tree-building using 1,000 ultra-fast bootstraps and automated substitution model optimisation with ModelFinder. Attempts to assess the functionality of genes based solely on genomic sequence can lead to the erroneous designation of functional genes as pseudogenes. Thus, for the purpose of this analysis, we consider genes as only the number of copies of a gene and do not speculate on functionality. Human, Xenopus and Drosophila BMP pathway genes are well-characterized and gene counts are not re-assessed here.

Gene family evolution was modelled using CAFE (v5.0.0)190. CAFE implements a stochastic birth and death model to estimate the number of gene gains and losses occurring at each node in a tree. The species set used for this analysis was identical to that used for phylogeny reconstruction with one exception: the nemertean Notospermus geniculatus was removed because the published genome has high levels of redundancy, which would result in inaccurate estimates. To run CAFE, an ultrametric version of the above species tree was calculated using pyr8s from iTaxoTools191 to implement r8s192. OrthoFinder was used as above for the orthology assignment. Principal component analysis was performed on the orthogroup-species matrix using R.

Macrosynteny analysis

Macrosynteny was compared between L. anatina and species representing four other phyla with relatively conserved genomic organizations: B. floridae (Florida lancelet, Chordata), P. maximus (scallop, Mollusca), L. longissimus (bootlace worm, Nemertea) and O. fusiformis (Annelida). Proteomes for these species were obtained from NCBI. SyntenyFinder (v1.0)85 was then used to run OrthoFinder to identify single-copy orthologues across the five genomes. To place the macrosynteny results in the evolutionary context of chromosome evolution across bilaterians, genes were assigned to the 24 bilaterian ancestral linkage groups (ALG)81 by homology. SyntenyFinder then runs the R package RIdeogram (v0.2.2) to visualise the outcome with idiogram plots and Oxford dot plots.

The genomic locations of BMP pathway genes, Wnt ligands, and Hox genes in each of the five species were identified with BLAST (v2.14.1+). We then determined the ALGs present on the chromosome hosting each gene using the above macrosynteny analysis. Conserved associations with ALGs were inferred from the genes’ chromosomal locations: if a gene was located on a chromosome containing the same ALG in all five species, then it was considered to have a conserved ALG association. To determine whether the rate of conservation of ALG associations is elevated in BMP pathway genes, Wnt ligands, and Hox genes compared to the background rate of all genes in the genome, the rate of conserved associations of developmental genes with a specific ALG was compared to that of a random sample of 100 single-copy orthologues using a chi-square test (Supplementary Data 37). The random sample was created by running OrthoFinder to identify single copy orthologues across the five species and then using the shuf command in bash with option -n 100 to choose 100 random rows of the output single copy orthogroups file. The presence or absence of conserved ALG associations for these genes was then determined by the same method as above.

Manipulation of BMP signals

BMP signalling manipulation experiments were deployed to reveal the function of the BMP pathway during L. anatina development. For the BMP(-) condition, two small molecule inhibitors, LDN-193189 (Stemgent 04-0074)86 and K02288 (Tocris 4986)87, were used to block the BMP pathway by inhibiting type I receptors for BMP ligand proteins. Two independent inhibitors were used to confirm that the observed phenotypes resulted from genuine BMP signal perturbation rather than being specific to a particular inhibitor. Though both molecules are highly specific to BMP receptors as opposed to TGF-β, and K02 has been reported to have no effect on TGF-β-induced Smad2 phosphorylation, they may have a weak off-target effect on Activin signalling87,193,194. To check the potential for off-target effects on TGF-β and particularly activin signalling in our system, we annotated the two activin/inhibin-like genes and other TGF-β-like genes in L. anatina (Supplementary Fig. 30) and quantified their expression during embryonic development. All TGF-β-like genes including the two activin/inhibin-like genes are either not expressed or expressed at negligible levels at the time of experimentation (Supplementary Fig. 31), suggesting that this pathway is not active at this stage of development and therefore that the effects of LDN-193189 and K02288 are not attributable to disturbed activin signalling.

Exogenous recombinant mouse BMP4 (mBMP4; R&D Systems 5020-BP) was used to overactivate the BMP pathway88. This condition is referred to as BMP(+). Manipulations were applied for two different durations, from the early blastula stage (5 h post fertilization, hpf) to either the late gastrula (10 hpf) or the one-pair-cirri larval stage (24–27 hpf). Two doses (high and low) of each manipulator were applied to the late gastrula experiments: LDN-193189 (4 and 8 μM), K02288 (250 and 500 nM) and mBMP4 (100 and 200 ng/mL). One dose was applied to the one-pair-cirri larva experiments: LDN-193189 (2 μM), K02288 (250 nM) and mBMP4 (200 ng/mL). Controls were run for each manipulator with the vehicle only (bovine serum albumin [BSA] or dimethyl sulfoxide).

Functional transcriptomics

RNA sequencing was used to explore the impacts of BMP signal manipulation experiments on the transcriptome. RNA was extracted for three biological replicates of each condition (45 samples, Supplementary Data 21) using TRIzol and sequenced using the Illumina HiSeq 4000 platform (total 855 million read pairs). After performing quality control with FastQC (v0.11.5) and trimming with Trimmomatic (v0.36), transcript abundances were quantified using kallisto (v0.43.0). Differential expression analysis was conducted using a bundled script in Trinity (v2.3.2)195. For datasets without biological replicates, edgeR was used with a fixed dispersion value of 0.1 to enable conservative inference196. For newly generated datasets with biological replicates, the voom method was applied to model the mean–variance relationship prior to linear modelling197. Transcripts were considered statistically differentially expressed with a false discovery rate (FDR) of less than 0.05. Pearson’s correlation coefficient was used to assess the similarity of each condition’s transcriptome. Putative BMP signalling downstream genes were identified based on their expression changes (at least fold-change > 2; P < 0.01). Genes that were upregulated in the BMP(+) condition and downregulated in the BMP(-) condition were classified as BMP-upregulated. Conversely, genes that were downregulated in the BMP(+) condition and upregulated in the BMP(-) condition were classified as BMP-downregulated. The transcriptomes of samples receiving high versus low doses of each manipulator were indistinguishable using hierarchical clustering (Fig. 4a) so are combined and treated as a single condition in subsequent analyses. Gene ontology analysis was performed on upregulated and downregulated gene sets using reciprocal best hits from BLAST searches against the Swiss-Prot database from UniProt198. To focus on the most robustly differentially expressed genes, we limited the dataset to genes that were differentially expressed in both the late gastrula and one-pair-cirri-larva experiments. Differentially expressed gene sets were searched for genes known to be involved in dorsal–ventral patterning in Xenopus.

Transcriptome age index analysis

The phylostratigraphic age of each gene in the L. anatina genome was first estimated using GenEra (v1.2.0)199,200. GenEra reduces biases in age assignment by using DIAMOND (v2.1.8) to search the entire NR database (Supplementary Data 38). To improve the resolution of age assignment in the case of L. anatina, we also added one genome assembly for each animal phylum that does not have a RefSeq annotated genome (Supplementary Data 39). This is especially important within the Spiralia, where key phyla like Bryozoa and Phoronida are unrepresented. Using the estimated gene ages, the R package myTAI (v1.0.1.9000)201 was used to calculate the transcriptomic age index122 for several L. anatina transcriptomic datasets. These were (i) a developmental time course from fertilised egg to two pair cirri larva (ii) L. anatina adult tissues (adductor muscle, lophophore, ovary, mantle and pedicle) and (iii) L. anatina late gastrula and larval stages with BMP signalling manipulation.

Bacterial-cloning-free riboprobe preparation

To prepare DNA templates for RNA probe synthesis, a bacterial-cloning-free protocol was developed to maintain both rapid preparation and target specificity (Supplementary Fig. 32). This protocol initially uses gene-specific primers (F1 and R1, designed from transcriptomes) to amplify target sequences in the first PCR. The amplified sequences are then ligated into the pGEM-T Easy Vector to attach an RNA polymerase promoter site. Inserts in the reverse direction (antisense) to the T7 promoter site were further amplified using T7 and gene-specific forward nested primers (F2) in a second PCR. This nested PCR ensures both the specificity of the PCR products and the correct transcriptional direction of the inserts. The second PCR products were used for sequencing validation and in vitro transcription. No cloning-based screening is required to select target clones, allowing the entire procedure to be completed within one day.

Total RNA from various embryonic stages was extracted using TRIzol and cleaned up with the RNeasy Micro Kit (Qiagen 74004). cDNA synthesis was performed using the SuperScript III First-Strand Synthesis System (Thermo Fisher Scientific 18080051). Target sequences were amplified in a first PCR using gene-specific primers (F1 and R1) and EmeraldAmp GT PCR Master Mix (TaKaRa RR320A). The PCR conditions were: 94°C for 2 min, followed by 30 cycles of 94°C for 15 seconds, 55°C for 15 seconds, and 72°C for 1 min 30 seconds, with a final extension at 72°C for 2 min. PCR products were analysed on a 1% agarose TAE gel stained with SYBR Safe DNA Gel Stain (Thermo Fisher Scientific S33102), and remaining products were purified using the QIAquick PCR Purification Kit (Qiagen 28104).

Purified PCR products were ligated into the pGEM-T Easy Vector (Promega A1360) using 5 μl of 2X Rapid Ligation Buffer, 0.5 μl of pGEM-T Easy Vector (25 ng), 3.5 μl of PCR product, and 1 μl of T4 DNA Ligase in a 10 μl reaction, incubated for 1 h at room temperature. A second PCR was performed using T7 and gene-specific forward nested primers (F2). The PCR conditions were: 94°C for 2 min, followed by 35 cycles of 94°C for 15 seconds, 55°C for 15 seconds, and 72°C for 1 min 30 seconds, with a final extension at 72°C for 2 min. PCR products were analysed on a 1% agarose TAE gel, and target bands were purified using the Wizard SV Gel and PCR Clean-Up System (Promega A9282). Primer sequences for probe synthesis are provided in Supplementary Data 40.

In vitro transcription was performed using T7 RNA polymerase (Promega RP2075) and DIG RNA Labeling Mix (Roche 11277073910). The reaction was incubated for 2–4 h at 37°C, followed by DNase I treatment for 15 min at 37°C. RNA probes were precipitated with 50 μl of nuclease-free water, 30 μl of 7.5 M LiCl, and 300 μl of 100% ethanol, chilled at -20°C for at least 30 min, and centrifuged at 16,000 × g for 20 min at 4°C. The RNA pellet was washed with 300 μl of cold 80% ethanol, air-dried, and resuspended in 25 μl of nuclease-free water. RNA probes were then diluted with 25 μl of 100% formamide and stored at -20°C or -80°C.

In situ hybridization for localising gene expression

Embryos were fixed overnight at 4°C with 4% paraformaldehyde (PFA; Electron Microscopy Sciences 15714) in filtered seawater. Post-fixation, embryos were washed with filtered seawater, dehydrated in 100% methanol, and stored at -20°C. For rehydration, embryos were transferred from methanol into baskets with 40-μm nylon mesh, immersed in PBST (0.1% Tween 20 in PBS) in a 24-well plate, and incubated for 10 min with 1 ml per well. Permeabilization was performed in PBSN (1% NP-40 and 1% SDS in PBS) for 10 min, followed by PBSTX (0.2% Triton X-100 in PBS) for 10 min. Optional bleaching was done in 2% H2O2 in PBST at room temperature for 30–60 min under direct light. Embryos were washed in PBST for 5 min, repeated three times, and rinsed with wash buffer (50% Formamide, 5X SSC, 1% SDS, 5 mM EDTA, and 0.1% Tween 20). Prehybridization was conducted in hybridization buffer (50% Formamide, 5% Dextran sulfate, 5X SSC, 1% SDS, 1X Denhardt’s, 100 μg/ml Torula RNA, 5 mM EDTA, 0.1% Tween 20, and 50 μg/ml Heparin) at 60°C for at least 1 h with slight rocking, with plates sealed in a plastic bag or covered with plastic wrap to prevent evaporation.

For hybridization, the hybridization buffer with probe (1:100–1:250, >100 ng/ml) was preheated at 70°C for 5 min, then applied to the embryos and incubated at 60°C overnight with slight rocking. Post-hybridization washes were done at 60°C for 15 min each, with three washes in wash buffer, followed by washes with wash buffer + 2X SSCT (1:1), 2X SSCT, 0.2X SSCT, and 0.1X SSCT. Embryos were then rinsed with wash buffer at room temperature for 15 min per wash. Blocking was performed in MAB blocking buffer (100 mM Maleic acid, 150 mM NaCl, 2% Blocking reagent [Roche 1096176], 10% Sheep serum, 0.1% Tween 20, and 0.2% Triton X-100) at room temperature for at least 1 h. Primary antibody incubation was carried out overnight at 4°C using anti-DIG-AP solution (Roche 11093274910) in MAB blocking buffer. Embryos were washed in MABTX (100 mM Maleic acid, 150 mM NaCl with 0.1% Tween 20 and 0.2% Triton X-100) five times for 20 min each, followed by two washes in TMN buffer (100 mM NaCl, 50 mM MgCl2, 100 mM Tris-Cl, 0.05% Tween 20) for 5 min each.

For the chromogenic reaction, embryos were transferred to clean wells to avoid staining from impurities. The reaction was carried out using BM Purple (Roche 1442074) in the dark without rocking, ranging from 20 min to several days depending on the probe. The reaction was stopped by washing in PBST, fixing in 4% PFA/PBST for 20–30 min, washing in 100% ethanol for 5–10 min, and rinsing in PBST. Finally, embryos were mounted in 70% glycerol in PBS with 0.1% NaN3 and stored overnight for full immersion at 4°C, protected from light. The embryos were imaged using a Zeiss Axio Imager Z1 or Axio Imager A1 upright microscope with DIC optics.

Immunostaining and imaging

Immunostaining was performed to identify the presence of phospho-Smad1/5 (pSmad1/5; BMP signalling readout) and phospho-Histone H3 (pHistone H3, indicator of mitotic cell at metaphase) during L. anatina development under control and BMP signal manipulation conditions. Embryos from early cleavage to larval stages were fixed in 4% PFA in filtered seawater, followed by dehydration in chilled methanol and storage at -20°C. For immunostaining, embryos were first rehydrated with PBST for 10 min and then subjected to permeabilization using PBSTX for 30 min. To block non-specific antigens, embryos were treated with 3% BSA in PBST for at least 1 h. They were subsequently incubated with either rabbit anti-phospho-Smad1/5 (1:200; Cell signalling 9511S) or rabbit anti-phospho-Histone H3 (Ser10) (1:100; Millipore 06-570) antibody in 3% BSA in PBST overnight at 4°C. Alexa Fluor goat anti-rabbit secondary antibody (1:400; Invitrogen A-11037) was used for signal visualization of the primary antibodies. For chitin detection, a fluorescein-conjugated chitin-binding probe (1:200; NEB P5211S) was applied. Nuclei were stained with Hoechst 33342 (1:1000 dilution from a 10 mg/mL solution; Invitrogen H1399), and cytoplasmic membranes were labelled using CellMask Deep Red (1:2000; Invitrogen C10046). Imaging was performed on a Zeiss LSM 780 or LSM 710 confocal microscope.

Statistics and reproducibility

Two-sided t tests were used for statistical comparisons as described in the figure legends. P < 0.05 was considered statistically significant. Exact P-values are provided in the figure legends or in the Supplementary Data. Whole-mount in situ hybridization and immunostaining experiments were performed using embryos from at least three independent fertilization batches with similar results. For each in situ hybridization or immunostaining condition, at least 30–50 embryos were examined and representative images are shown. For transcriptomic analyses, three biologically independent samples per condition were used. No statistical method was used to predetermine sample size. Embryos from the same fertilization batch were randomly allocated into experimental groups. Blinding was not applicable, as BMP signalling manipulation produced visible phenotypic effects on embryonic development.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "Lewin et al. present a chromosome-scale genome of the brachiopod Lingula anatina and show that 4 out of 10 Lingula chromosomes are results of fusion-with-mixing. By comparing the genomes of a mollusc, an early-branching annelid and Lingula they find evidence for lineage-specific fusions. Lewin et al. then focus on BMP signaling components, since the role of BMP signaling in the DV patterning of Spiralia is insufficiently well-understood. They show that BMP2/4 and Chordin are expressed on the opposing sides of the DV axis (dorsally and ventrally, respectively), and demonstrate a pSMAD1/5 gradient with a dorsal maximum, similar to the situation in the mollusc Ilyanassa and non-spiralian bilaterians. Given the diversity of the ways spiralians pattern their DV axis, this manuscript provides important support for the hypothesis that Chordin-BMP-dependent DV patterning was an ancestral feature of spiralian development shared not only by Spiralia but by Bilateria as a whole. It is an important paper, and I strongly recommend publication, however, there are several points, which need to be resolved first.\n\nMajor points\n1. Both, in the Introduction and in the Discussion, the Authors give a misleading, oversimplified view of the BMP signaling pathway, which is reflected in their Fig. 2b and Fig. 2f, and which causes some obvious misinterpretation of the results. The Authors seem to think that BMP2/4 (=Dpp) is the only BMP responsible for the DV patterning. This is not the case. First, BMPs function as dimers, and usually the DV-patterning-relevant ligands are the BMP2/4/BMP5-8 heterodimers. So, BMP5-8 is not an \"alternative BMP paralog\" relevant only for the leech. It is involved in DV patterning across Bilateria and in the directive axis patterning in anthozoan Cnidaria. All the models Lewin et al. mention in the Introduction use BMP5-8 together with BMP2/4. In Xenopus, BMP7 (belongs to the BMP5-8 group) is expressed ventrally together with BMP4; in Drosophila, Screw (=BMP5-8) is expressed ubiquitously, while Dpp is dorsal; in leech, BMP5-8 is dorsal and BMP2/4 is ubiquitous etc.\nSecond, the Authors repeatedly call ADMP a BMP antagonist, but it is not. It is a BMP ligand (Reversade and De Robertis, 2005). And this misunderstanding probably made Lewin et al. overlook their own most striking result – they showed a “BMP seesaw” in a spiralian without noticing it. In anthozoan Cnidaria as well as in ambulacrarian and chordate deuterostomes, the second body axis is patterned by a \"BMP seesaw\". Both ends of the second body axis express BMPs, but BMP signaling takes place on the end of the axis, which is opposite to Chordin. For example, in Nematostella, one end of the directive axis expresses BMP2/4, BMP5-8 (BMPs) and Chordin, and the other end expresses GDF5-like (BMP). BMP signaling is on the GDF5-like side. In sea urchin, BMP2/4, ADMP1 (BMPs) and Chordin are expressed ventrally, and ADMP2 (BMP) - dorsally. BMP signaling is dorsal. In the frog, BMP4 and BMP7 (BMPs) are ventral, BMP2, ADMP (BMPs) and Chordin are dorsal. BMP signaling is ventral. Lewin et al. see exactly the same in a spiralian, with dorsally expressed BMP2/4 and BMP5-8 (BMPs), and ventrally expressed ADMP, BMP3 (BMPs) and Chordin and pSMAD1/5 opposite to Chordin. This is absolutely fantastic and has to become the focus of the revised manuscript. In the revision, the Authors must show pSMAD1/5 stainings and in situs of Chordin and all BMPs expressed at gastrula (BMP2/4, BMP5-8, BMP3, ADMP) in controls and upon up- and downregulation of BMP signaling on one main text figure – not scattered between the main text and extended figures – and make a really big evolutionary point about the “BMP seesaw” being potentially the ancestral DV patterning mechanism in the Discussion (stressing the similarity to the mollusc data doi:10.1093/molbev/msab322 and to non-spiralian data). In contrast, Bambi and CV2 in situs are, of course, good to have – as examples of conserved BMP modulators expressed in the “high BMP signaling” region across species, but they can easily go into the supplement. However, in situ quality needs to be improved, especially for BMP2/4, BMP5-8 and CV2, where signal-to-noise ratio is unacceptable. The Authors should maybe consider clearing their embryos with one of the recently published protocols e.g. https://www.science.org/doi/10.1126/sciadv.aba0365 to improve the quality both of the in situs and of the morphological description shown in Extended Data Fig. 7.\n2. If the Authors publish a chromosome-scale genome, they should also make it easy to use for the community. I understand that this is a lot of work, but a publicly accessible Lingula genome browser with a built-in BLAST option rather than a link to a bioproject with raw reads would be great.\n\nOther comments\n1. Line 81-88. The \"typical spiralian\" argument is not ideal. Among high-level spiralian taxa, spiral cleavage exists in mollusks and annelids, however, as far as I am aware, no lophophorates (bryozoans, phoronids etc) exhibit spiral cleavage. Flatworms or nemertines also do not have spiral cleavage. Most gnathiferans (chaetognathes, rotifers etc.) also do not have spiral cleavage, but gnathostomulids have it. Given that, I do not think we can state with certainty that spiral cleavage was lost multiple times or gained multiple times in Spiralia (as the Authors do by saying “possibly reflecting a reversion to the ancestral state”). On the other hand, some spirally cleaving molluscs use BMP for DV patterning, and spirally cleaving annelids maybe use it as well - for example Owenia has Chordin and BMP2/4 and likely uses it for the DV patterning. In general, I suggest that instead of explaining why other models are not ideal, the Authors should argue that with the increased sampling across the spiralian phylogeny, we get a better understanding of the ancestral mode of BMP signaling in Spiralia.\n2. Fig. 1c. How did the Authors get the list of BMP signaling component orthologs for Fig. 1c? Was it BLAST or reciprocal BLAST or something else? How certain are the Authors that lophotrochozoan BMP3 is an ortholog of chordate BMP3? I may have easily missed something, but if I am not mistaken, the only published non-chordate putative BMP3 comes from the ctenophore Mnemiopsis, however, its phylogenetic position as a BMP is not robust (compare Fig. 2 in doi: 10.1371/journal.pone.0024152 and Supplementary Fig. 2 in doi:10.1242/dev.141507). If BMP3 is indeed always present in Lophotrochozoa, it is an important piece of information. A phylogenetic tree with BMP ligand distribution across Bilateria (including ADMP, GDF5-6 etc.) would be a very useful supplementary figure. It would also be very useful to add Xenopus to Fig. 1c either in addition or instead of human, since this is the animal authors constantly refer to in the text.\n3. Fig. 1d. I assume the Authors did not do this for space reasons, but they should add Owenia and Linneus data on Fig. 1d. Having them separately on Extended Fig. 6a-b is very confusing, especially since the order of the chromosomes is always different for each species on each figure. If having amphioxus, Pecten, Lingula, Owenia and Linneus on one panel is impossible (although I do not see why it should not be!), and the Authors want to keep their Pecten + Lingula + “1 more species in the middle” comparisons, they should fix the order of Pecten and Lingula chromosomes, adjust the third species to this fixed order and put all these figures as panels on the main text figure.\n4. For ALG association of BMP genes, the result is quite striking, but I may have overlooked a careful explanation of how random genes were selected. It is important and must be clearly explained in the main text. Also, the Authors should comment on whether the chromosomal linkages of BMP signaling components are preserved in highly scrambled genomes (Drosophila, C. elegans, clitellate annelids, tunicates etc.) or genomes which underwent more translocations than, say amphioxus or Pecten - e.g. human genome, Streblospio benedicti, Carcinoscorpius rotundicauda.\n5. Lines 239-240 and Fig 2a. The Authors write that “The co-option of BMP–Chordin signalling pathway into dorsal–ventral axis patterning is a cornerstone of the bilaterian body plan”. While I fully agree that BMPs and Chordin are important, I am not sure about the co-option. Co-opted from doing what? Co-option means that BMP and Chordin were doing something else in the bilaterian ancestors and then were re-deployed in DV patterning. If their function in DV was re-deployment, what was the ancestral function? Also in the Fig. 2a, the Authors suggest that “The BMP-Chordin axis predates the origin of bilaterians and, therefore, the dorsal–ventral (DV) axis.” I am not sure how the Authors make this conclusion. What do they call a BMP-Chordin axis? A body axis patterned by BMP and Chordin? Their tree on Fig. 2a suggests that it evolved prior to the split of Cnidaria and Bilateria. In that case, cnidarian second body axis (directive axis of anthozoans) is homologous to the DV axis of Bilateria, and it should be therefore called DV axis as well. However, it is not entirely certain whether directive and DV are homologous or evolved convergently. If it is convergence, then Lewin et al. cannot say that \"BMP-Chordin axis predates the origin of Bilateria\" - in this case, these are different body axes, which just happened to use BMP signaling for patterning. I suggest removing the term “co-option” and carefully rephrasing the legend.\n6. Fig. 2b. I understand that the Authors tried to keep things simple, however, not only do BMPs signal as dimers (see major point 1) but also the transduction of the signal is performed not by a type I BMP receptor but by a heterotetrameric complex containing two type I BMP receptors and two type II BMP receptors. Additionally, BMPR-I is located within the cell in the schematic, however the receptor complex should be positioned on the cell surface instead.\n7. Fig 2c. What is the Author´s definition of the “main BMP signalling ligands, mediators, and modulators”. I would recommend either including the full heatmap of L. anatina BMP pathway components shown in extended data Fig. 8a or restricting the listed components in Fig 2c to only the BMP ligands (including ADMP) and Chordin, as these are the focus of the paper.\n8. Fig. 2d-e. It is not entirely clear why these have to be main text figures. It is nice to show that Lingula SMAD1/5 is really a SMAD1/5 (Fig. 2d), and it is excellent that the Authors demonstrate that the epitope recognized by the anti-pSMAD1/5 antibody used in the paper is present in Lingula SMAD1/5 (Fig. 2e), however, this are “technical details”, which should go into the supplement. Also, the Authors start using the term R-Smad without explaining what it is. They must write that Smad1/5 is the R-Smad in the BMP signaling (and Smad2/3 is the R-Smad in the Nodal/TGFb signaling) - otherwise unclear.\n9. Fig. 2f. I find this very confusing. First, the signal is normally transduced not by a Smad1 homotrimer as the Authors depict on the image but by a heterotrimer containing two Smad1/5 molecules and one Smad4 (Co-Smad) because such heterotrimer is energetically more stable than a homotrimer. Second, I do not understand why they point out the L3 loop - they never discuss this in the text, and this is not the site where their antibody binds. L3 should either be removed from the image and legend or elaborated upon in the text. If no additional discussion of the Fig. 2f content is added, I suggest its removal from the paper as it does not contribute significantly to the main narrative.\n10. Fig. 3a. Why were two treatment durations used? This is not clear and must be explained in the main text. In addition, it must be made clear in the main text or figure legend of each experiment throughout the entire manuscript which treatment regime was used.\n11. Fig 3, Fig 5, extended data figure 9. Sample size and phenotype penetrance should be provided for all immunostainings and RNA in situ hybridisations (i.e. X embryos out of a sample of Y develop the phenotype shown on the representative image).\n12. Fig 3c-f. It is mentioned briefly in the main text that manipulation of BMP signalling delayed gastrulation. Do the authors have data regarding the extent to which gastrulation is delayed? Is this the case for all individuals? Do the animals eventually complete gastrulation?\n13. Fig. 3g-j. Why are 2 different doses of inhibitor/recombinant BMP used? Which dosage were the animals shown in panels g-j exposed to? It must also be made clear in the figure legend or main text whether the animals used in experiments shown in fig. 3c-f and 3k-n were treated with high or low doses. This information must also be made clear for all experiments throughout the manuscript that made use of the BMP inhibitors/recombinant BMP.\n14. Fig. 3k-n. This needs better (larger and better annotated) pictures. Also, the Authors had two types of treatment - do they see the chitin phenotype in both or only in the treatment spanning longer? They have to comment on whether there are any differences between the treatments.\n15. Lines 376-377. GO terms are useful but also very imprecise. I think that calling these genes “neuronal” is a BIG stretch. I do not know which alk the Authors mean - for me Alk2, Alk3/6, Alk4/5/7 are different type I BMP or TGFb receptors, which cannot be called \"neuronal\". FoxB is primarily an endomesodermal marker. Lhx1/5 has functions in the neurons but also in the development of the lymphoid tissue. Netrin is an axon guidance molecule but has other functions as well. Nkx2.4 is expressed in the brain, but predominantly in the testes. Other genes the Authors mention here can probably be called neural without too many additional disclaimers. I do not doubt that in some context all these genes are expressed in neurons, but I think that this bit needs some very careful re-phrasing. E.g. \"dominated by genes demonstrating either exclusively neuronal expression or known to be expressed in neurons in specific developmental contexts in other model organisms\". I understand that in situ validation of neuronal expression of all these genes in Lingula maybe too much for this paper (although it would be great to have it), but without it the argument has to be tuned down and carefully phrased.\n16. Fig 4g. BMPR should be positioned on the plasma membrane, not inside the cell.\n17. Fig. 5a. The Authors should replace BMP5 with BMP7 in the Xenopus ventral centre, add BMP2, and re-color ADMP in the dorsal centre (ADMP is not a BMP antagonist! It is a BMP!).\n18. Lines 480-490. The section on TAI analysis is confusing. It is unclear how the results of this analysis indicate that BMP signals are crucial for gastrulation (Line 488). Is it because of the GO terms associated with the gastrulation process? The delayed gastrulation observed in Fig. 3 is more convincing evidence that BMP signalling is important for gastrulation. In addition, the description of the role of BMP signalling as “crucial” for gastrulation is too strong, as it is previously stated on line 319 that BMP manipulations only delay gastrulation, and not that BMP inhibition prevents it.\n19. Supplementary fig 5a. In the text, late blastula and mid gastrula stages are highlighted as stages with the highest TAI score but the Authors refer to a figure highlighting a different stage – the late gastrula, which has the lowest TAI.\n20. Supplementary figure 5b. It should be made clear what “BMP related” genes are. Are these genes that are members of the BMP family, or are they components of the BMP pathway?\n21. Lines 491-499. This needs to be re-written once my Major point 1 has been addressed. Also more data needs to be presented to confirm that BMP2/4 and BMP5-8 expression on the dorsal side is suppressed by BMP signaling. It is counterintuitive, and if true, very interesting. However, it is not sufficiently well demonstrated. High quality in situs and QPCR are necessary upon up- and downregulation of BMP signaling.\n22. Lines 512-513. Hemichordates and sea urchins are not \"basal deuterostomes\". Together, they are Ambulacraria and a sister group to Chordata. If the Authors wish to take sea urchins and hemichordates separately, then Chordata would be an a earlier branch of Deuterostomia than either Hemichordata or Echinodermata.\n23. Fig.6. The BMP gradient depicted (pink) should be much narrower. The width of the pink shaded region better reflects the region of Dpp expression than the nuclear pMad gradient (Raftery and Sutherland, 2003).\n24. Fig 6. In the legend, the blue region on each embryo cartoon is described as the “anti-BMP domain”. This is inaccurate/vague since this area expresses its own BMP ligands and it should just be referred to as the chordin/sog expression domain.", "reviewer_2": "BMP signaling is involved with establishing the dorsal-ventral axis and a region of neuroectoderm in several taxa in deuterostomes and ecdysozoans. However, experiments testing the function of BMP signaling in these processes across spiralian taxa have found differing levels of involvement. To understand whether BMP signaling was ancestrally involved in forming the dorsal-ventral axis including the neuroectodermal domain in bilaterians, the authors examined BMP function in the spiralian Lingula anatina, a brachiopod that develops by radial rather than the ancestral spiral cleavage. Unlike some other spiralian taxa, L. anatina appears to have some of the same regulatory logic of the BMP signaling gradient in chordates, including an anti-neural function for BMP. The authors also generated a chromosome-level assembly for L. anatina and reconstructed the complement of BMP pathway components for multiple spiralian taxa. The paper is well-written and appropriate for a broad audience, the data are of high quality, and the experiments and findings are important for understanding how the dorsal-ventral axis including the neuroectoderm evolved in bilaterians.\n\nMain comments/suggestions\nThe authors argue that the radial cleavage present in L. anatina may be a ‘reversion to the ancestral state’. However, this still raises the question of the evolutionary history of BMP signaling in spiralians. Since spiral cleavage was likely ancestral for the clade containing brachiopods, annelids, and mollusks, then BMP signaling was first used in the context of spiral cleavage. It would be interesting for the authors to comment on how they think BMP signaling was used to form the dorsal-ventral axis in spirally-cleaving embryos and was then re-evolved to form the dorsal-ventral axis in a radially-cleaving embryo. Why do you think there is such variability in the involvement of BMP signaling and dorsal-ventral axis formation across different spiralian taxa? I think the discussion points on lines 563–570 are very interesting, and begin to get at this this question, but I would like the authors to comment on this in the context of gain and loss of spiral cleavage. In L. anatina, BMP signaling does not appear to be active (pSMAD1/5 levels) during cleavage stages, unlike other taxa with spiral cleavage, but instead is active in the blastula. To really demonstrate a lack of BMP signaling during cleavage, I think it’s important to test function of BMP during this time window (see additional comments for Fig. 2 below).\n\nThe regulation of bambi, cv2, admp, and chordin by BMP signaling in L. anatina (Fig. 5) is interesting. The authors suggest that this is indicative of a BMP signaling gradient controlling fates along the dorsal-ventral axis as in chordates, and that this system may have been present in the last common ancestor of the two groups (e.g., lines 450–451). An alternate interpretation of the functional transcriptomic data is that a gradient of BMP signaling was ancestral and that regulatory logic is maintained in L. anatina, irrespective of what developmental program (e.g., fates along the dorsal-ventral axis) this signaling system controlled. I think that ventralization or dorsalization (or loss/gain) of multiple fates after drug treatments or BMP protein treatment in L. anatina needs to be demonstrated more clearly to support the hypothesis of an ancestral function in dorsal-ventral axis formation. The change in expression domains of bambi, cv2, admp, and chordin after perturbing BMP signaling does not necessarily indicate a shift in fates along the forming dorsal-ventral axis, especially since many of the other Spemann-Mangold organizer genes were not found to be regulated by BMP signaling in the same direction as in Xenopus (e.g., lines 476–478).\n\nIn L. anatina, does the neuroectoderm form on one side of the trunk ectoderm during gastrulation, and are the differentially-regulated genes expressed in the neuroectoderm? It looks like the references for the sentence on lines 375–377 are for other animals. For Fig. 4, I would like to see spatial expression by in situ hybridization of the some of the differentially expressed neural genes identified (e.g., elav, pax3/6, soxB1) in controls and before and after manipulation of BMP signaling.\n\nAdditional comments and questions\nIntro\nline 52 Capitalize “ADMP” since it’s an abbreviation for “anti-dorsalizing morphogenetic protein”.\nlines 55–56 “Basal deuterostomes” would be the nodes before the terminal/extant taxa. I think you mean “early-branching deuterostomes”.\n\nResults and Figures\nFig. 2\nIn panels m, n, and o, there are some cells where the pSMAD1/5 labeling appears to only be on one side of the nucleus. Is this common for pSMAD1/5 labeling and if so, what does it indicate? I’m asking about the cells with a circular domain of labeling, not the cells with punctate spots of pSMAD1/5.\nThe view in panel k is not labeled (av?).\nIs vegetal is down in panel l?\nFor panel o, I am not certain what I am looking at. Why is there a stripe of pSMAD on the ventral side of the larva? Is this a lateral view similar to Fig. 3k? A diagram of the tissues including the forming mantle and pSMAD in an early larva would be helpful.\nWhy do you think there is no pSMAD1/5 at cleavage stages even though there is maternal bmp5-8 and smad1/5 transcript and ubiquitous bmp receptor expression? Did you try blocking BMP signaling with the two drugs during cleavage stages (i.e., a drug pulse during early cleavage and assessment of resulting larval phenotypes)? Similarly, if you add BMP protein during cleavage stages, do you see pSMAD during these stages? Is it possible that there is a lower level of BMP signaling during cleavage that is not detected by the cross-reactive antibody? Is bmp5-8 transcript spatially restricted to a subset of cells during cleavage?\n\nFig. 3\nFor the BMP manipulation experiments, the long time-window went from 5 hpf to ~24-27 hpf. What was the evidence that the BMP protein or the BMP inhibitors were active for that whole time in L. anatina? Was the BMP protein or BMP signaling inhibitors only applied once at 5 hpf? It would be useful to show pSMAD1/5 at early larval and 1 PCL stages for the long treatment times.\nFor the two BMP protein doses (low and high) and the two drug doses, did you test whether these doses resulted in different levels of pSMAD1/5 in L. anatina? Are the data shown in Fig. 3d–n after high dose treatment for the protein and the drugs? I don’t see data comparing the different doses except effect on pH3, and the change in number of pH3+ cells between the low and high treatments does not look very dramatic compared to the shift from no treatment to the low treatments.\nHow was the overall morphology including the dorsal-ventral axis affected by up or downregulating BMP signaling? You state that gastrulation was delayed and that mantle lobe folding was disrupted, but a description of additional phenotypes would be useful for understanding how manipulating BMP signaling affected development and formation of fates along the dorsal-ventral axis.\nWhat is the view in Fig 3n?\nIt would be helpful to include a sentence explaining where in the BMP signaling pathway LDN-193189 (LDN) and K02288 (K02) act. Since both drugs inhibit multiple ALK receptors, what is the evidence that the dose you are using is selective for BMP signaling and not Activin signaling in L. anatina? Did you test the effects on pSMAD1/2? Similarly, for the BMP protein treatments, did you test whether you were cross-activating Activin signaling (SMAD1/2)? It’s important to demonstrate the selectivity for perturbation of only BMP signaling since Activin signaling has been shown to be an organizer signal in at least one other spiralian.\n\nFig. 4\nFor the functional transcriptomic experiments, did you find differences in which genes were up and down-regulated between the 100 and 200 ng/mL BMP treatments or between the low and high drug concentrations? Do you think these low and high concentrations are both maximally upregulating/blocking BMP signaling, or did you find evidence for BMP acting as a morphogen in L. anatina?\n\nFig. 5\nIn panel e, the z-scores for some genes (e.g., shh, twsg1, bmp5, bmp4, fst…) seem to vary across replicates within a treatment. Is this due to low expression levels for these genes?\nIn Extended Fig. 9, the up-regulation of bmp2/4 expression in the ectoderm after LDN treatment is not very clear to me in the animal shown.\n\nDiscussion\nBe careful of inferring gene function in L. anatina based on what a gene homolog does in chordates. For example, in lines 515–516, you discuss “ventralising factors, particularly chordin”. However, the function of chordin in L. anatina has not been tested to my knowledge.\n\nMethods\nline 589 It looks like multiple tissues were collected for genome sequencing, but only mantle, lophophore, and adductor muscle were used. Was the DNA in the other tissues degraded?\nline 596 For genome sequencing, did you assay DNA integrity using BioAnalyzer in addition to assessing purity with a NanoDrop?\nline 601 & 738 Did you do any column purification after RNA extraction by TRIzol?\nline 742 Why was edgeR used versus DESeq2, and what were the criteria for using a dispersion parameter of 0.1?", "reviewer_3": "I co-reviewed this manuscript with one of the reviewers who provided the listed reports. This is part of the Nature Communications initiative to facilitate training in peer review and to provide appropriate recognition for Early Career Researchers who co-review manuscripts.", "reviewer_4": "The role of BMP signaling in patterning the dorsoventral body axis has been described in a diversity of species, including deuterostomes, ecdysozoans and spiralians, and is widely considered to be the ancestral bilaterian state. While initial studies in a handful of model organisms, primarily in vertebrates and insects, suggested BMPs are required to inhibit neural induction in early embryogenesis, recent work in emerging model systems, particularly those found in the Spiralia, suggest a more complex story where neural induction is independent of BMP signaling in multiple species (e.g. Lambert et al., 2016; Webster et al., 2021; Webster & Meyer, 2024).\n\nIn the manuscript “Brachiopod genome unveils the evolution of the BMP-Chordin network in bilaterian body patterning” Lewin et al., present a chromosomal-scale genome for the brachiopod Lingula anatina. This genome improves on a previous scaffold-level assembly for the same species and extends observations from this previous publication from a number of these authors that BMP patterning may be involved in shell formation in brachiopods. While the chromosome-level genome will prove to be a key resource for future work in this clade, the study largely focuses on the role of BMP signaling in this species. The authors employ the genome to describe the complement of BMP signaling pathway members in this species as well as other spiralians. The experiments presented, including pharmacological approaches and RNA seq, confirm a role for BMP in patterning the dorsoventral body axis, which has been suggested in the previous genome description and was previously described in other spiralians, including two other brachiopod species, T. transversa and N. anomala (Martin-Duran et al., 2016).\n\nBased on differential expression in RNAseq datasets generated from pharmacologically up and down regulating BMP signaling, the authors suggest that BMP also inhibits neural development in L. anatina. They find dozens of genes to be differentially regulated upon manipulation of BMP signaling, including a number of neural-associated genes. Based on this finding, the authors conclude that BMP is required to inhibit neural development, a role well described in vertebrates and arthropods. However, the experimental evidence for this assertion is not sufficient, and indeed, expression patterns for several neural genes contradicts this hypothesis. For example, elav is unchanged or even slightly downregulated in BMP(-), while soxB1 is upregulated in both BMP(-) and BMP(+) embryos. Furthermore, the (differential) expression of neural genes should also be confirmed using in situ hybridization assays to visualize where these genes are expressed in relation to the BMP signal and potential changes in expression domain after treatment.\n\nThroughout the manuscript, the authors discuss the roles of BMP in DV-axis patterning and in neural development as one common feature, with one being dependent on the other. However, recent work in other spiralians indicates that these can be two distinct roles that can be present independently of each other. The authors imply that because they see an effect on dorsoventral patterning, a role in neural development must be the ancestral state and vice versa. This discussion would benefit from clear delineation and test of these different roles, and would greatly benefit from an attempt to differentiate these functions experimentally. Indeed, recent work in spiralians suggests a more complex story for the role of BMPs in neural development, with functional studies in both annelids and molluscs demonstrating neural development independent of BMP signaling. This needs to be discussed in more detail, taking into account different scenarios like convergent evolution or loss of traits.\n\nSome additional concerns and comments:\n1. The role of Chordin: The title and abstract of the manuscript refer to the role of the BMP-Chordin network. While the manipulation of BMP signaling seems to affect expression levels of Chordin (whether directly through transcriptional regulation or indirectly through changes in size of the ventral domain is also not addressed), the actual role of Chordin in body axis patterning and neural development is not directly investigated.\n2. Expression of BMPs and pSMAD signal\nUsing antibody staining, the authors detect phosphorylated SMAD1/5, a read-out for BMP signaling activity, at the dorsal side of the embryo (Figure 2). In Extended Data Figure 9, they also show the expression of bmp2/4, bmp3 and bmp5-8, using ISH. This expression data should be included in Figure 2. The expression data for bmp2/4 and bmp5-8 expression raises concerns. While there is a dorsal signal that likely overlaps with pSMAD1/5 signal, bmp2/4 and bmp5-8 also seem to be expressed more globally throughout the embryo. Can the authors provide an explanation for this, or ISH images that are easier to interpret for non-brachiopod specialists?\n3. Role of BMP in cell proliferation\nIn Figure 3, the authors describe an effect of BMP manipulation on cell proliferation. This observation should be discussed in more detail in the manuscript. What cell types are affected by this role? How does this compare to known roles of BMP signaling in other systems? How could this be related to other observed phenotypes?\n4. Role of BMP in shell formation\nIn Figure 3, the authors describe a role of BMP signaling in establishing the embryonic shell fields. They conclude that BMP is required to divide the shell field into two valves, as the fields seem to be fused in BMP(-) embryos. However, in BMP(+) embryos the shell field seems to be expanded and experience over-folding, rather than lost as might be expected if BMP indeed inhibits shell field formation. How can these results be explained?\n5. RNA sequencing\nThe RNA sequencing data require further confirmation and discussion: 1) How much overlap is there between DEGs in BMP(+) and BMP(-) embryos? Are the same pathways/genes affected? 2) Is there evidence in the RNA sequencing data for the described role in cell proliferation and shell formation/folding? What are affected signaling pathways? 3) How reliable is GO term analysis in this species? 4) As mentioned above, the changes in neural gene expression require further confirmation/investigation (see above). 5) The authors describe two different incubation windows used for BMP manipulation. It is not clear which one of these was used for RNA sequencing.\n6. The Lophophorata hypothesis: The authors present a phylogenetic analysis in Fig1b that suggests that brachiopods are in a clade with phoronids and bryozoans to form the “Lophophorata”. While this is a clade proposed in historic analyses, recent phylogenomic work suggests that bryozoans may be sister to entoprocts or cycliophorids, which the authors suggest may be the result of including orthologs with lower phylogenetic signal. However, these results are difficult to interpret as the authors did not include entoprocts or cycliophorids in their analysis. Could the authors search for their gene set in the deposited data for these species to provide support for their hypothesis? Or, could they leverage their new chromosome-level assembly could be to search for shared derived syntenic relationships between brachiopods and bryozoans that might support their hypothesis?\n7. Terminology\nThe authors make references to “basal deuterostomes” throughout the manuscript. This should be corrected to “basally-branching deuterostomes” or referred to as the phylogenetic group specifically (e.g. ambulacrarians).\n8. Does the proposed karyotype comport with expectations in this group?\n\nOther minor comments:\nMain Figures\nFigure 1, Panel D: Why is the Hox cluster pointed out specifically? It is hardly mentioned in the text. The lower half of this panel should be discussed and explained in the figure legend.\nFigure 1, Panel E: This panel requires more explanation.\nFigure 2, Panel D-F: These could go into the supplement.\nFigure 2, Panel G-O: Labeling of body axis and potentially an embryo schematic would be helpful here. Expression of BMPs should be added to this Figure.\nFigure 3, Panel K-N: Overlaying the chitin staining with the pSmad staining would be helpful to confirm that BMP signaling is active in the Chitin-negative regions.\nFigure 3, legend c-f: What do the empty arrowheads indicate?\nFigure 3, legend k-n: What does the arrowhead indicate?\nFigure 4, Panel A-H: Panels A and G-H require additional explanation or should be presented in a more intuitive way. Panels B-D do not provide crucial information and could be moved to Supplement.\nFigure 4, Panel A: Pearson correlation between ctrl and BMP(-) conditions show surprisingly little difference. How can this be explained?\nFigure 4, Panel A: What are the different samples that each column and row represents?\nFigure 4, Panel C-D: This could go into supplementary figures.\nFigure 5, Panel B: Body axes should be labeled for easier understanding. For admp and bambi the difference between control and BMP(-) embryos is not very clear. A quantification or at least a better representative image should be provided.\nExtended Data\nFigure 7: This figure is not mentioned in the text.\nFigure 9, Panel A: ISH images are difficult to interpret. While a bmp2/4 signal is present in the dorsal part of the embryo where pSmad was detected in Figure 1, it also seems to be expressed in other parts of the embryo?\nFigure 9, Panel B: Where is bmp2/4 marked in the schematic?\nFigure 9, Panel C: This requires further investigation. Can this be confirmed with the RNAseq data? Is the effect on chordin expression direct or indirect because the ventral domain is lost?\nFigure 9, Panel D: What experimental evidence is the distinction between ectoderm and endomesoderm based on?\nFigure 9: Please provide more explanation for the arrows and pathways in the Figure legend.\nSupplementary Data\nFigure 4: Why was the tree only produced for the Smad family? At the very least this should be produced for all major proteins discussed in this study (i.e. BMP family, Chordin, selected neural genes, dorsoventral markers). Also: add description of species abbreviations.\nText\nPage 1, Line 30: “We uncover a BMP signaling gradient…” - The presence of a signaling gradient was never actually discussed/shown. For BMP only mRNA expression was shown, and not the distribution of protein. We appreciate the challenges of working with emerging model systems, but this should be clarified in the text. The pSmad antibody could be used as a read-out for a signaling gradient, but this would need to be discussed in more detail.\nPage 2, Line 52: “...In the fruit fly Drosophila, the BMP gradients is completely inverted…” - This is true for invertebrates in general.\nPage 9, Line 244: “...Extended Data Fig. 7…” - This should probably refer to Extended Data Fig. 8?\nPage 11, Line 306: “...underscoring a deeply conserved BMP gradient…” - See above. The use of the word gradient should be justified.\nPage 11, Line 310: “...we applied [...] inhibitors and [...] proteins for two durations…” - The reasoning for selection of these stages/durations should be explained.\nPage 11, Line 314: “...LDN-193189 (LDN) and K02288 (K02)...” - What is the difference between these two inhibitors? Why were two used?\nPage 11, Line 316: “...led to a loss of dorsal pSmad1/5…” – From the images presented, it looks like the dorsal signal persists, and additional signal is gained throughout the rest of the embryo.\nPage 11, Line 319: “...altered gastrula morphology…” - The authors should describe how gastrula morphology is altered and what can be concluded for the role of BMP signaling in the process.\nPage 13, Line 371: “Gene ontology (GO) analysis [...] highlighted the ‘nervous system development’...” - How reliable is GO term analysis in this species? Phylogenetic trees should be presented, at least for the main protein families of interest to confirm this.\nPage 17, line 477: “...organiser genes cer1, noggin and shh are BMP-upregulated when they might be expected to be downregulated…” - The authors should provide confirmation of this with in situ staining. This also merits further discussion.\nPage 17, Line 491: “Finally, we investigated how BMP ligand expression related to the BMP signaling gradient…” - This should be shown in Figure 1 or 2 and not the extended data."}, "author_response": "Point-by-point response to the reviewers’ comments\nWe would like to begin by thanking all four reviewers for their detailed and thoughtful\ncomments on the initial version of our manuscript. We have addressed nearly all of the\npoints raised and have substantially revised many parts of the work as a result. In particular,\nkey improvements include:\n1. creation of a Lingula genome browser with integrated BLAST search capabilities;\n2. inclusion of a more complete and detailed description of BMP signalling in the\nIntroduction;\n3. reworking of the narrative to focus on the identification of a BMP seesaw in a\nspiralian;\n4. clarification of the distinct roles of BMP signalling in dorsal–ventral patterning versus\nneural induction;\n5. repetition of in situ hybridisation experiments using an optimised protocol, yielding\nclearer results;\n6. addition of new in situ experiments revealing neural gene responses to BMP\nmanipulation, supporting transcriptomic data;\n7. construction of gene trees for all gene families characterised in this study.\nBelow, we provide a point-by-point response to each comment detailing the changes made."}, {"round": 2, "reviewer_comments": {"reviewer_1": "Lewin et al. have significantly improved the manuscript, and I am happy the changes. There are some textual things, which, in my opinion, should be addressed before the acceptance, but I do not think another revision round is necessary. The data are excellent, and once the authors fix the text a little, the paper should be accepted.\nHere are some comments in the order of appearance in the text.\n1. Throughout the manuscript, the authors state that both in vertebrates and invertebrates, the neurogenic domain forms on the “low BMP signaling” side of the DV axis and mention the ancestrality of the centralized nervous system in Bilateria as a proven fact. “CNS in low BMP zone” also marked on the root of their tree on Fig. 6. This is a misleading oversimplification. Insects with their ventral CNS do not represent the variety of NS types in ecdysozoans, and not every spiralian has a neat ventral nerve cord like Platynereis. Also ambulacrarian deuterostomes do not have a CNS. There are many invertebrate groups with much less centralized – sometimes even diffuse – nervous systems, but a dorsal BMP signaling maximum is clearly an ancestral feature of Bilateria. This automatically means that some neural tissue in such “non-CNS” animals must be forming in the “high BMP” area. See https://www.sciencedirect.com/science/article/pii/S001216061730475X and https://www.biorxiv.org/content/10.1101/2025.06.08.658475v1.full It is very important that the authors discuss the negative regulation of the potentially pro-neural markers in Lingula, however, they also must put this into context. Lewin et al. also need to elaborate on what is known about the embryonic origin of neurons in Lingula and the extent of centralization of its nervous system. Does the SoxB2-expressing area of the Lingula gastrula really represent neurectoderm? Or do Lewin et al. think that it is neurectoderm because it expresses SoxB1, SoxB2 and Pax6? How centralized is the NS of the Lingula larva and adult?\n2. Line 89-90: The authors talk about the “variability of ligand localisation” – please replace with “variability of ligand expression”.\n3. In the new version of the manuscript, the authors replaced “basal deuterostomes” with “early-branching deuterostomes”. Ambulacrarians are as early-branching as the chordates. I would just stick to “ambulacrarian deuterostomes”.\n4. Line 463: Lewin et al. write about “neurulation in spiralians”. Spiralians do not neurulate. Neurulation is a chordate-specific term.\n5. In their discussion of the organiser (Lines 622 and further), the authors raise a question of whether the organiser may have been present outside chordates. The blastopore lip organiser is clearly much older than chordates – it exists in cnidarians. See https://www.nature.com/articles/ncomms11694 However, molecularly, cnidarian organiser resembles the Nieuwkoop center more than the Spemann-Mangold organiser. If Lingula has a blastopore lip organiser, it would be important to know whether only Chordin-expressing side is capable of axis induction.", "reviewer_2": "I appreciate the efforts the authors made to revise the manuscript; my concerns and questions have been adequately addressed. The revised focus on BMP signaling---description of a BMP ‘seesaw’, separation of d-v axis formation and neural induction as potentially separate events, and a more nuanced description of BMP signaling and putative evolutionary scenarios throughout----has greatly improved the manuscript. I also think the addition of new gene expression data for ‘neural’ genes strengthens the hypothesis that BMP signaling affects neural fate in L. anatina. I also appreciate the clarifications throughout the manuscript in terms of methodology, reporting of results, and interpretation of those results. I think this manuscript makes important contributions to our understanding of the function of BMP signaling during development in Spiralia and evolution of BMP signaling across Bilateria and should be published.\nI have a few minor questions, but I don't need a response to these before publication.\nIs there a chordin-like homolog in addition to chordin L. anatina?\nAre there distinct elav1 and elav2 homologs in L. anatina?\nI am still not convinced that the neural tissue is increased/decreased after perturbation of BMP signaling in L. anatina. The expression of elav after manipulating BMP signaling is not very strong.\nFig. 2e – How do you know which side is ‘dorsal’ from a lateral view (is the future dorsal side morphologically distinct from the ventral side)?"}, "author_response": "We would like to take this opportunity to once again thank the reviewer for their time and effort in reviewing our manuscript. Their comments were instrumental in many improvements we were able to make to the manuscript, and we are delighted that they are satisfied with the result."}]

Decision letter