Sequence alignments were performed on amino acid sequences deduced using MUSCLE 3.8.425 (Edgar, 2004) within the Geneious Prime suite of software (https://www.geneious.com). Candidate homologous sequences for Marchantia genes were obtained from Phytozome (https://phytozome-next.jgi.doe.gov/), Sol Genomics Network (https://solgenomics.net/) or NCBI (https://www.ncbi.nlm.nih.gov/), and are provided in Supplementary Data S1-1 . Maximum likelihood phylogenetic trees were generated on the alignments using PhyML 2.2.4 (Guindon and Gascuel, 2003) with default settings.

For plant transformation, M. polymorpha L. spores were obtained as described by (Albert et al., 2018). Plant lines were maintained asexually through the propagation of gemmae, either plated directly on 0.5× Gamborg’s B5 medium (Duchefa Biochemie, Haarlem, the Netherlands; 1% (w/v) sucrose, 1% (w/v) agar) or onto a sterile filter paper disc covering the medium. Standard culture conditions were 25°C, 16 h photoperiod and 30 µmol m^-2 s^-1 light intensity provided by cool white, fluorescent tubes. For carbon and nutrient deprivation experiments, gemmae (four independent transgenic lines and >50 biological replicates per treatment) were plated onto water agar medium solidified using 1% agar and grown for 10 days either in the light/dark or complete darkness.

Three candidate Marchantia genes related to SL pathways were targeted: two genes for DWARF 27, MpD27–1 (Mp6g03970) and MpD27–2 (Mp6g01750), and a CAROTENOID CLEAVAGE DIOXYGENASE 7 (MpCCD7/Mp2g03280). A gene encoding a protein known to regulate gemmae biology, SUPPRESSOR OF MORE AXILLARY GROWTH2-LIKE (MpSMXL/Mp3g06310), was also targeted. CRISPR/Cas9-based genome editing was performed as previously described (Albert et al., 2018). Agrobacterium-mediated transformation was used for transforming spores with a construct containing four guide RNAs per gene (guide sequences are provided in Supplementary Data S1-2 ). The four 20 bp guides/gene (N17VVR: protospacer adjacent motif) (Doench et al., 2014) were designed using Geneious Prime and synthesized by GenScript Biotech (NJ, USA) as a polycistronic sequence, consisting of 20 bp guide sequences fused to a modified sgRNA scaffold (Dang et al., 2015) separated by glycyl-tRNA sequences. The synthesized sequences were subcloned into pMpU6ENTR and recombined into the binary vector pMpGE010 (Kawamura et al., 2022) using LR Clonase II enzyme mix (Thermo Fisher Scientific). Following transformation, PCR amplification and DNA sequencing were used to find mutations within the target genes ( Supplementary Figure S1 , Supplementary Data S1-3 ). At least five independent mutant events were shown for each of the targeted genes. Mutants were found in the initial transgenics (T0) because the main life-stage of Marchantia is gametophytic and, thus, haploid. Mutant lines were propagated through gemmae (G1 generation) and re-sequenced to ensure establishment of non-chimeric lines, as gemmae are derived from single cells. Plants used for experiments were G2 generation or later. For most experiments, four independent lines were used, and all analysis were performed using the parental line as a control.

To assay for variation in gemmae numbers, we assigned numbers to the gemma cups on the thallus of 8-week-old wild-type plant. The Marchantia thallus grows outwards from the meristem at the thallus branch tip, showing new gemma cups as it grows. Thus, the oldest cups are at the thallus base and the youngest next to the meristem. Gemma cups at the thallus base were called Gemma Cup 1 (GC1) and successive cups were numbered outwards to the tip of the thallus branch (GC7). For comparisons between transgenic lines and wild-type, GC3 and GC4 were used. To quantify gemmae numbers, gemma cups were removed from the thallus and counted using a magnifying lens.

To quantify gemmae growth during the standard day/night growing conditions, we photographed the growing gemmae and then measured gemmae area by using Image J (Schneider et al., 2012). The germination of gemmae held in the dark was scored by counting the number of elongated gemmae and gemmae that developed rhizoids.

The wild-type control and three transgenic lines Mpd27–1#16, Mpd27–2#9 and Mpccd7#17 were grown on half-strength Gamborg’s B5 medium for 5 weeks. For uniform tissue sampling, a 6 mm diameter disc encapsulating GC3 and GC4 was harvested, weighing approximately 160 mg per biological replicate (with three biological replicates per genotype). Total RNA was extracted using the Spectrum™ Plant Total RNA Kit (Sigma Aldrich) following the manufacturer’s instructions. RNA-sequencing library construction and sequencing was performed by Beijing Genomics Institute (Shenzhen, China) using DNBSEQ Eukaryotic Strand-specific mRNA library preparation and the DNBseq platform. The resultant reads were mapped onto M. polymorpha subsp. ruderalis genome version 6.1 (https://marchantia.info/) using HISAT2 (Kim et al., 2019) with default parameters. Read counts were then calculated by Feature Counts (Liao et al., 2014). Feature counts generated were further filtered for genes with lower counts and replicates where genes are not expressed at all using edgeR (filterByExp, with default settings, ensuring removal of genes with lower than 10 counts in all samples). This left a total number of 13,408 genes across all samples. The featurecounts were also processed for K means clustering (top most variable 2000 genes, 6 clusters) and Hierarchical clustering analysis using the iDEP application through the subtraction of each gene’s mean expression level (Ge et al., 2018).

The filtered featurecounts were then processed for pair-wise differential expression (DE) analysis using edgeR (QLFTest) (Robinson et al., 2010). Files for each comparison were then filtered for genes with FDR <0.05 and logFC +-1.5. These files are presented for each comparison and an Upsetplot was performed on all these treatments (https://upsetplot.readthedocs.io/en/stable/formats.html). Raw counts data are described in Supplementary Data S1-4 , while differentially expressed genes (DEGs) are summarized in Supplementary Data S1-5 .

Ten milligrams of powdered freeze-dried thallus tissue from 8 weeks old Marchantia plants was mixed with 1 mL of acetone, homogenized using a vortex mixer for 30 s, and shaken at 4°C for 30 min to allow extraction of the lipophilic compounds. After centrifugation at 10,000 ×g for 10 min, the supernatant was collected. The extraction was repeated on the pellet using the same extraction solvent. Both supernatants were combined, evaporated until dry under a gentle stream of nitrogen, and resuspended in 1 mL of ethanol. The resulting extract was filtered through a 0.45 µm syringe filter and stored at −20°C prior to UPLC analysis. The quantification was carried out using a Waters Acquity UPLC system (Milford, MA, USA) equipped with a Photodiode Array Detector and a Mass Single-quadrupole Detector (QDa, Waters). An aliquot of 2 μL was injected onto an Acquity UPLC BEH C18 column (2.1 × 100 mm, 1.7 μm particle size, Waters) at 40°C, with a flow rate of 0.5 mL min^-1. The solvents were: Solvent A: 0.1% formic acid in water, Solvent B: 0.1% formic acid in acetonitrile, and Solvent C: 0.1% formic acid in isopropanol (all v/v). The gradient elution program was: 0 min: 30% A, 66% B, 4% C, 0.5 min: 30% A, 66% B, 4% C, 10 min: 18% A, 72% B, 10% C, 13 min: 2% A, 78% B, 20% C, 17 min: 2% A, 78% B, 20% C, 17.1 min: 30% A, 66% B, 4% C, 20 min: 2% A, 78% B, 20% C. Carotenoids were detected at 420 nm and quantified as lutein or β-carotene equivalents using a six-point calibration curve and lutein and β-carotene calibration standards prepared from commercially available sources. Data acquisition with the QDa mass detector was undertaken in positive mode using the following conditions: total ion current (TIC) between m/z 100–1000 Da, capillary voltage: 1 kV, probe temperature: 600°C, and cone voltage: 5 V. Compounds were confirmed in selected ion recording (SIR) mode.

Droplet Digital PCR (ddPCR) was performed to validate the DEGs obtained from the RNA sequencing experiment. We used the BIO-RAD droplet digital PCR system to quantify transcripts of the target genes, with MpACTIN serving as the reference gene. The primer sequences are provided in Supplementary Data S1-3 .

Gene homologs for two steps of the SL biosynthetic pathway have been identified from Marchantia, DWARF27–1 and DWARF27–2, (MpD27–1/Mp6g03970 and MpD27–2/Mp6g01750) and CAROTENOID CLEAVAGE DIOXYGENASE 7 (MpCCD7/Mp2g03280). However, no gene sequences for two pivotal enzymes responsible for catalysing the final stages of this pathway, MORE AXILLARY GROWTH 1 (MAX1) and CAROTENOID CLEAVAGE DIOXYGENASE 8 (CCD8), have been identified in previous work (Delaux et al., 2012; Walker et al., 2019; Kodama et al., 2022) ( Figure 1B ). To examine the relationship of the MpD27–1, MpD27–2 and MpCCD7 candidate genes to those of other species, we conducted a phylogenetic analysis of D27 and CCD7 deduced amino acid sequences from representative species of liverworts, lycophytes, gymnosperms, and flowering plants ( Figure 1C ). Using iterative BLAST searches, a single CCD7 homolog was identified in Marchantia (MpCCD7), along with four homologous D27 genes (MpD27–1, -2, -3 and -4). Of these, MpD27–1 and MpD27–2 were chosen for further analysis because they belong to the same clade. The phylogenetic analysis supported the identification of a single CCD7 ortholog and four closely related D27 homologs in Marchantia ( Figure 1C ).

Although the function of these genes in M. paleacea has not been characterized yet, it possesses a functional SL pathway (Kodama et al., 2022). Therefore, based on the close sequence identity, MpD27–1, MpD27–2 and MpCCD7 were hypothesized to perform functions related to SL metabolism and chosen for further functional analysis. Based on homology analysis these genes seem to regulate biosynthesis of the SL precursor, 9-cis-10′-apo-beta-carotenal. Considering that ancestors of liverworts had SLs (Kodama et al., 2022), we are referring to these genes in Marchantia as SL precursor biosynthesis genes, despite the absence of detectable SL production in Marchantia. This categorization is based on sequence similarity, probable evolutionary relatedness, and likely functional conservation.

To investigate the biological functions of MpD27–1, MpD27–2 and MpCCD7, CRISPR/Cas9 Marchantia knockout mutants were generated for each gene, using Agrobacterium-mediated transformation of spores. The resultant mutations ranged from single nucleotide insertions or deletions to larger modifications, with deletions extending up to 370 bp and insertions reaching 231 bp ( Supplementary Figure S1 ). Most of the genetic changes were predicted to result in truncated proteins. Marchantia is dioicous, so spores for transformation must be generated from a cross of male and female lines, which can generate heterogeneous genetic background among sporelings (Ishizaki et al., 2016). Thus, to isolate the effects of the mutations from genetic background variability, multiple knockout lines were identified for each target gene. Four distinct mutant lines were consequently employed as biological replicates for all experiments (except transcriptomic analysis), in addition to the parental line. Additionally, based on initial analysis of the mutant lines suggesting a phenotypic change in gemmae biology and gemmae gemination on the mother plant, a CRISPR/Cas9 knockout mutant for SUPPRESSOR OF MORE AXILLARY GROWTH2-LIKE (MpSMXL/Mp3g06310) was generated, as this gene regulates karrikin signaling and gemmae numbers in Marchantia (Komatsu et al., 2023) and provided a useful comparative control.

In tissue culture conditions, the gemmae of wild-type plants typically remain within the gemma cups. However, Mpd27–1, Mpd27–2, and Mpccd7 mutants displayed a higher percentage of gemma cups with gemmae outside of the cups compared with the wild-type control ( Figures 2 , 3A ). Gemmae, upon dispersal from the cup, initiate growth, a phenomenon also observable in SL mutants, as evidenced in Supplementary Figure S2A . Previously, it was reported that mutants of MpSMXL exhibited an increase in gemma number, retarded thallus growth and suppressed gemma dormancy in the dark (Mizuno et al., 2021; Komatsu et al., 2023). Given that karrikins are thought to potentially mimic the bioactivity of SLs, we hypothesized that an increase in the number of gemmae per cup could be a driving factor behind this gemmae escape phenomenon. The Mpsmxl#11 mutant did not exhibit the visible increase in gemmae outside of the cup observed in the SL precursor mutants ( Figures 2 , 3A ). However, to appropriately test the hypothesis, we quantified the gemmae number for both the SL precursor mutants and wild-type. To facilitate this analysis, we first examined how the spatial position of the cups influenced gemmae numbers in wild-type genotype.

A numbering system was used based on the relative age of gemma cups. The oldest cup, situated at the thallus base, was designated as GC1 (Gemma Cup1) and younger cups closer to the apical notch were numbered consecutively from GC1 (reaching GC7), as illustrated in Figure 3B . No significant differences were observed in the number of gemmae within cups positioned at various locations along the thallus ( Supplementary Figure S2B ). For consistency, subsequent analyses concentrated on GC3 and GC4, with the gemmae from both cups combined for quantification. The Mpd27–1, Mpd27–2, and Mpccd7 mutants had no significant difference in gemmae numbers from wild-type ( Figure 3C ). This suggested that the observed gemma cup release phenotype in the SL precursor mutants was not due to an increase in gemmae count. Thus, MpD27–1, MpD27–2 and MpCCD7 may have a role in suppressing release of gemmae from the gemma cup, thereby reinforcing gemmae dormancy while on the parent plant.

Given that MpD27–1, MpD27–2 and MpCCD7 influence gemmae dispersal, a process that typically initiates gemmae germination under favorable conditions, we investigated whether the genes also impact gemmae germination and growth. We compared gemmae growth under nutrient- and carbon-starved conditions either in the presence or absence of light. The mutant lines grown in the light did not show any significant differences in gemmae growth compared with wild type, as determined by comparative surface area ( Supplementary Figure S3 ).

The emergence of rhizoids is an indicator of the release of gemmae dormancy: gemma without rhizoids are considered dormant and those with visible rhizoids categorized as non-dormant (Eklund et al., 2015, 2018). Additionally, gemma germination and growth in the dark is typified by elongation of one side of the gemma (Mizuno et al., 2021). To evaluate germination of dark-grown gemmae, elongation and/or the presence of rhizoids were measured in mutants and the wild-type line. This assessment was impractical for light-grown gemmae, as all initially displayed gemmae growth and rhizoid emergence, until development was halted due to nutrient scarcity. Consequently, area-based metrics were exclusively used to gauge growth of light-grown gemmae. In dark conditions, all four lines of Mpd27–1, Mpd27–2 and Mpccd7 exhibited significantly higher numbers of gemmae with elongation compared with wild-type line ( Figure 4 ). Additionally, the majority of the Mpd27–2 mutant lines and all the Mpccd7 mutant lines exhibited an increase in rhizoid production ( Supplementary Figure S4 ). For Mpsmxl#11 gemmae, growth appeared like wild-type line. The experimental evidence shows that MpD27–1, MpD27–2 and MpCCD27 are required to keep gemmae dormant under unfavorable (nutrient depleted and dark) conditions.

To find changes in mRNA abundance associated with gemma release, germination, and growth, comparative transcriptome analysis of Mpd27–1, Mpd27–2 and Mpccd7 mutants and wild-type control was conducted, focusing on evaluating the thallus tissue surrounding GC3 and GC4 ( Figure 5A ). Principal Component Analysis (PCA) demonstrated that biological replicates within each genotype were grouped closely ( Figure 5B ), whereas distinct groupings were observed among different genotypes, affirming the reliability of the data and the distinctiveness of phenotype between mutant and wild type. The raw counts and results of DEG analysis are given in Supplementary Data S1-4, S1-5 .

Hierarchical clustering of the transcriptomic data grouped the top 2,000 most variable genes into two main clusters ( Supplementary Figure S5 , Supplementary Data S1-6 ). Cluster 1 included the transcript data from Mpccd7 and Mpd27–2, while cluster 2 was composed of wild-type and Mpd27–1. This indicated that Mpd27–2 and Mpccd7 exhibited similar gene expression patterns, whereas Mpd27–1 differed from other genotypes. Using K-means clustering yielded the expression patterns into notable six clusters ( Figure 5C , Supplementary Data S1-7 ). Notably, clusters 4 and 2 exhibited the most contrasting profiles based on genotype-based differences. Furthermore, cluster 5 and 6 consisted of genes associated with the genotype’s sex.

To identify the genes that were differentially expressed in a pair-wise fashion between SL mutants and wild type we kept genes that showed log2FC of ± 1.5 at FDR <0.05. Mpd27–1 had 421 downregulated and 117 upregulated genes, Mpd27–2 had 160 downregulated and 45 upregulated genes, and MpCCD7 showed 139 downregulated and 37 upregulated genes ( Figure 5D , Supplementary Data S1-7-S1-9 ).

Employing UpSet plot on this dataset ( Figure 5E , Supplementary Data S1-8 ), identified that 59 genes were differentially regulated across all three mutants. Specifically, 27 genes showed differential expression between Mpccd7 and Mpd27–2, while 37 genes differed between Mpd27–1 and Mpccd7. Among the detected DEGs, with the FDR/fold-change criteria used of <0.05/>1.5-fold, 388 out of 538 were associated only with Mpd27–1. In contrast, far fewer DEGs were found only in the Mpd27–2 and Mpccd7 datasets with these criteria, being 65 out of 205 and 53 out of 176, respectively.

To enhance the clarity of our findings, we concentrated on analyzing DEGs that exhibited consistent changes across all three mutants. Specifically, we focused on genes predicted/known to have biological functions crucial to essential plant developmental processes including photosynthesis, stress adaptation, and terpenoid metabolism ( Figure 6 , Supplementary Data S1-9 ). Notably, genes encoding early light induced (ELI) proteins (including Mp4g18580 and Mp4g18590), peroxidases (e.g., Mp5g01640 and Mp3g07970), an α/β hydrolase (including Mp2g15700), and late embryogenesis abundant (LEA) proteins (including Mp4g14840 and Mp4g14870) were predominantly downregulated. Conversely, Mp7g00860/MpERF10, an ethylene response factor, demonstrated increased transcript levels. Genes associated with terpenoid metabolism exhibited variable responses. For example, Mp3g13150/MpTPS6 was significantly upregulated with a log2FC of 3.3, whereas Mp3g21760/MpTPSL18 was downregulated with a log2FC of -5.9 ( Figure 6A , Supplementary Data S1-9 ). Interestingly, the transcript abundance of Mp4g19090/MpTPSL23 was higher in Mpd27–2 and Mpccd7 mutants compared to the Mpd27–1 mutant ( Supplementary Data S1-9 ).

The DEGs associated with photosynthesis, hormones, and stress adaptation were shared between Mpd27–1 and Mpd27–2 mutant RNAseq datasets ( Figure 6B , Supplementary Data S1-9 ). For example, genes including CAB/ELI (e.g., Mp4g18610 and Mp4g18620), GERMIN-LIKE (e.g., Mp5g00880, Mp5g00920, and Mp5g13870), and ABA-responsive e.g., DEHYDRIN (Mp6g15610) and LEA (e.g., Mp7g06630 and Mp4g14880) were downregulated. At the same time BETA-CAROTENE DIOXYGENASE (Mp4g01260) was upregulated in both mutants. Similarly, DEGs between Mpd27–1 and Mpccd7 were mostly associated with stress responses. For example, Mp4g02840 (induced by phosphate), Mp8g18620, (A-type lectin), and Mp5g05460 (known as fasciclin-like) were downregulated, whereas Mp6g18420 (a cytochrome P450), showed an increase in expression ( Figure 6C , Supplementary Data S1-9 ). The downregulated genes common between Mpd27–2 and Mpccd7 included CAB/ELI (e.g., Mp4g18560), Stigma specific (e.g., Mp5g11400) and PALMITOYL GLYCEROLIPID (e.g., Mp6g04480) ( Figure 6D , Supplementary Data S1-9 ). Conversely, genes such as MpERF17/Mp4g22280, TRIHELIX16/Mp3g10500), MpTPSL23/Mp4g19090 and LRAT/Mp3g10470 were upregulated. Taken together, these findings indicate that SL precursor mutants affected many genes regulating essential biological functions such as photosynthesis, ethylene signaling, stress responses, terpenoid and carotenoid metabolism.

The transcriptome profile unique to Mpd27–1 included many genes linked to stress response and cellular transport ( Figure 7A , Supplementary Data S1-9 ). The ten lectin-related genes, involved in defense and communication responses against biotic and abiotic stress, were suppressed, exhibiting log2FC from -1.6 to -6.4. For example, A-type lectins Mp1g18730 with a log2FC of -2.59 and B-type lectins Mp8g18670 with a log2FC of -6.46 were downregulated. Similarly, transport-related genes involved in transport of iron (including Mp2g25340 with a log2FC of -3.3), phosphate (including Mp4g16570 with a log2FC of -3.46), and nitrate (including Mp4g03070 with a log2FC of -5.6) showed decreased transcript abundance. In contrast, some phosphate transporters like Mp4g11070 were upregulated with a log2FC of 3.3. Furthermore, terpenoid synthase (Mp6g04630/MpMTPSL10) increased by a log2FC of 2, and LOX2 (Mp2g12180) was upregulated with a log2FC of 3.46.

In Mpd27–2, MpCCD7 (Mp2g03280) was upregulated 4-fold, suggesting that the mutant might be responding to a defect in a metabolite biosynthesis process regulated by MpD27–2 ( Figure 7B , Supplementary Data S1-9 ). This change could represent a compensatory mechanism for reduced flux within the pathway or be a result of feedback regulation. These findings suggest a potential link between MpD27–2 and MpCCD7 in the same metabolic pathway.

The Mpccd7 mutant showed 3.6-fold increase in transcript abundance of 9-CIS-EPOXYCAROTENOID DIOXYGENASE (Mp2g07800/MpNCED) ( Figure 7C , Supplementary Data S1-9 ). Other genes with notable reduced expression included SPORE COAT POLYSACCHARIDE (Mp3g15460) with a log2FC of -1.5 and FASCICLIN-LIKE ARABINOGALACTAN (Mp5g05390) with a log2FC of -1.7. In contrast, MALATE SYNTHASE (Mp8g16290), MpCYP704-like8 (Mp3g10820), and MpNBS-LRR5 (Mp3g09150) were upregulated.

Digital Droplet PCR (ddPCR) was used to validate the RNA-Sequencing findings, by examining six genes identified as altered in the mutant lines ( Figure 7D ). The ddPCR results confirmed the differential expression of the genes compared with the wild-type control. For example, MpCCD7 (Mp2g03280) expression was increased >2-fold in Mpd27–2, while MpERF10 (Mp7G00860) exhibited a more than 4-fold increase in Mpd27–2 and an over 8-fold increase in Mpccd7. In contrast MpCAB (Mp4G18580) was -5-fold downregulated in the Mpd27–1, -3.6-fold in Mpd27–2 and -2.45-fold in Mpccd7 mutants. These results are consistent with the RNAseq data, thereby validating the observed differential gene expression.

Given their orthologous relationship with Arabidopsis thaliana’s D27, MpD27–1 and MpD27–2, are predicted to function as β-carotene isomerases. This is further supported by the observed increase in transcript abundance of carotenoid associated genes in transcriptome analysis of SL precursor mutants. To investigate whether the loss of gene function impacts carotenoid composition, a targeted UHPLC analysis was conducted on wild-type and mutant genotypes ( Figure 8A ). Additionally, the observed decrease in transcript levels for photosynthesis-related genes, as shown in Figure 6A , prompted the quantification of chlorophyll content in the mutants ( Figure 8B ).

The carotenoids violaxanthin, all-trans- and cis-neoxanthin, lutein, α-carotene and β-carotene have been identified in Marchantia, and several unidentified carotenoids are also present (Takemura et al., 2014). The identification and quantification of carotenoids in the lines studied here showed a significant increase in the levels of cis-neoxanthin and lutein across all three SL mutants. Additionally, violaxanthin levels were higher in Mpccd7, α-carotene was elevated in Mpd27–2 and Mpsmxl, and β-carotene saw an increase in Mpd27–1. An unidentified carotenoid designated as “compound 4”, also exhibited a significant rise in Mpd27–1 and Mpd27–2. Similarly, the chlorophyll data indicated that SL mutants had increased levels of chlorophyll a, b and pheophytin in comparison to the wild-type control ( Figure 8B ).

Bryophytes, which encompass non-vascular plants like liverworts, mosses, and hornworts, can utilize both sexual reproductive organs, such as spores, and vegetative structures like gemmae and some part of the thallus for their propagation (Wyatt, 1982; Chopra and Bhatla, 1990; Buck et al., 2003; Whitaker and Edwards, 2010). Gemmae are small disc-like asexual propagules produced in specialized structures called gemma cups, present on the dorsal side of the thallus (Barnes and Land, 1907, 1908). Raindrops have the potential to expel gemmae from gemma cups (sometimes referred to as ‘splash cups’), allowing plants to grow and establish themselves away from the parental plant (Edwards et al., 2019). The gemmae contained within the gemma cup remain dormant and do not undergo maturation or growth until they are either expelled from the cup or until the maternal plant reaches the end of its life cycle. This suggests a restraint on the liberation and germination of gemmae from the maternal plant. While there is a growing body of information about the genetic factors controlling growth and germination of gemmae (Eklund et al., 2015; Eklund et al., 2018), the genes that regulate the release of gemmae from the cup in liverworts remain elusive.

Here we have shown that defects in functions of MpD27–1, MpD27–2 and MpCCD7 regulate gemma cup release ( Figures 2 , 3A ). These Marchantia mutants, when compared with wild type, displayed elevated levels of transcripts for ETHYLENE REPONSE FACTORS (ERF) orthologs, e.g., Mp7g00860/MpERF10 and Mp4g22280/MpERF17 ( Figures 6A, D ). The expression of Mp7g00860/MpERF10 was also validated through ddPCR analysis, providing clear evidence of increased MpERF10 levels in all three mutants ( Figure 7D ). Furthermore, the MBEX expression tool indicated higher expression of MpERF17 in gemmae cups (Kawamura et al., 2022). ERFs, belonging to the family of AP2/ERF transcription factors, regulate many processes in flowering plants, including seed germination and pod shattering (Pirrello et al., 2006; Chung et al., 2010; Müller and Munné-Bosch, 2015; Chandler and Werr, 2020; Hu et al., 2022; Liu et al., 2022; Ali et al., 2023). If we consider gemmae release to be analogous to seed release, we suggest that the SL biosynthesis genes in Marchantia could act to restrain the function of ERFs, thereby inhibiting the release of gemmae.

Mutants of SL precursors exhibited a notable decrease in the transcript abundance of EARLY LIGHT INDUCED (ELI) genes ( Figure 6A , Supplementary Data S1-9 ). Specifically, Mp4g18580 and Mp4g18590 showed down regulation across all three mutants, while Mp4g18620 and Mp4g18610 were down regulated in both Mpd27–1 and Mpd27–2, and Mp4g18560 experienced a decrease in Mpd27–2 and Mpccd7. Further independent support for the role of these genes in gemma biology came from gene expression profiles. Results obtained using the MBEX expression tool revealed that all four ELI/CAB genes showed increased expression in the gemma cup (Kawamura et al., 2022). These proteins, while not directly binding chlorophyll (Green et al., 1991), are categorized as part of the CAB family and recognized as early light-induced proteins (ELIPs). ELIPs are triggered by light stress and implicated in pigment biosynthesis and thylakoid membrane assembly, essential for protecting chlorophyll-protein complexes from light damage.

SLs impact the process of senescence by regulating various genetic components involved in photosynthesis (Mayzlish-Gati et al., 2010; Mashiguchi et al., 2021; Li et al., 2022; Schiphorst et al., 2022). For example, the tomato SL biosynthesis Sl-ORT1 mutant showed reduced chlorophyll content and lower levels of light-harvesting genes (e.g., CAB and RUBISCO) (Mayzlish-Gati et al., 2010). Similarly, the Arabidopsis SL signaling max2 mutant exhibited delayed senescence and reduced expression of ELIP1, ELIP2, and HY5 (Shen et al., 2007). The decreased expression of ELI genes in mutants of SL precursors, along with their connection to gemma and gemma cup development, suggests a role for ELI in SL-mediated gemmae development. Furthermore, mutants of SL precursors exhibited higher chlorophyll content compared to the control plants, suggesting that SL regulates senescence in Marchantia ( Figure 8B ). However, it is important to note that the increase in chlorophyll content in SL precursor mutants may be attributed to an increase in the number of growing gemmae.

Our data show that MpD27-1, MpD27-2, and MpCCD7 inhibit the germination and growth of gemmae under dark and nutrient-starved conditions ( Figure 4 , Supplementary Figure S4 ).

This is consistent with other studies on bryophytes that report enhanced tissues growth in mutants defective in the SL and karrikin pathways. For example, SL-deficient (Ppccd8) and SL-insensitive (Ppkai2Lgjm) moss mutants showed increased extension of moss filaments in the dark (Lopez-Obando et al., 2021). Similarly, karrikin signaling mutants, Mpkai2a and Mpmax2, showed enhanced growth of gemmae under dark conditions (Mizuno et al., 2021).

SL and ABA are known to interact with each other in regulating plant stress responses. For instance, rice mutants defective in SL production (e.g., max4/dwarf10 and max3/dwarf17) and SL perception (e.g., dwarf3) showed higher concentrations of ABA (Haider et al., 2018). Similarly, ABA-deficient tomato mutants exhibited reduced concentrations of SL (López-Ráez et al., 2008, 2010). Furthermore, the external application of the SL analog GR24 in Lotus japonicus led to an increased transcript abundance of ABA biosynthesis and ABA responsive genes (Liu et al., 2015).

In Marchantia, ABA contributes to the induction of gemmae dormancy. For example, ABA-insensitive Marchantia mutants were unable to establish and maintain gemmae dormancy (Eklund et al., 2018). The reduced expression of ABA-biosynthesis genes, e.g., 9-CIS-EPOXYCAROTENOID DIOXYGENASE (MpNCED/Mp2g07800), and ABA-responsive genes, including DEHYDRIN (MpDHN3/Mp6g15610) and LATE EMBRYOGENESIS ABUNDANT (LEA/Mp4g14840, Mp4g14870, Mp7g06630, and Mp4g14880) ( Figure 7 ) in mutants of SL precursors suggests a linkage between SL signaling, ABA production, and the control of gemmae dormancy and germination processes. This conclusion is further supported by evidence from the study by Kawamura et al. (2022), where gene expression profiles generated using the MBEX expression tool indicated that all four MpLEA genes were expressed in the gemma cup and sporeling. Additionally, MpDHN3 and MpNCED exhibited higher expression levels in sporelings. Our analysis revealed that MpNCED expression was notably decreased in the Mpccd7 mutant compared to the Mpd27–1 and -2 mutants ( Supplementary Data S1-5 ). This decrease correlates with enhanced rhizoid production observed in the Mpccd7 mutant lines ( Supplementary Figure S4 ).

Mpd7–1 showed reduced transcript abundance of many LECTIN genes, specifically A-type and B-type lectins ( Figure 7A ). Notably, gene expression profiles generated using the MBEX expression tool indicated that Mp8g18650 was highly expressed in gemmae cups (Kawamura et al., 2022). Lectins play a pivotal in establishing a symbiotic relationship between plants and different type of organism under various stress conditions (De Hoff et al., 2009; Jain et al., 2022). For example, under nitrogen stress, legumes secrete lectins from their roots to facilitate bacterial attachment, leading to nodule formation (Downie, 2010). Similarly, certain lectin genes, specifically induced by mycorrhizal associations, are categorized as AM-induced lectin-like genes (Frenzel et al., 2005, 2006). In lichens, an algal-binding lectin has been found to play a key role in establishing a symbiosis between the mycobiont (fungal component) and the photobiont (photosynthetic component) (Nazem-Bokaee et al., 2021).

SLs are crucial plant hormones that regulate plant architecture and act as signaling molecules exuded into the rhizosphere to foster symbiotic relationships, specifically with AM fungi. Although various species of Marchantia, including M. pappeana and M. paleacea, interact with AM fungi, M. polymorpha stands out as it does not form association with AM fungi. Despite this M. polymorpha possesses the capability to interact with a wide range of microbial communities (Poveda, 2020). For example, M. polymorpha was able to establish interactions with various pathogenic and beneficial endophytic fungi (Nelson et al., 2018; Nelson and Shaw, 2019). Similarly, M. polymorpha and M. paleacea interacted with diverse bacterial genera known for plant-growth promotion, exudate degradation, nitrogen fixation and disease-suppression (Van Damme et al., 2008; Alcaraz et al., 2018).

The observed decrease in lectin genes expression in the SL precursor mutant Mpd27–1 suggests that MpD27–1 might be involved in establishing interactions with microorganisms through lectins during biotic and abiotic stress responses.

In SL precursor biosynthesis mutants, we observed an increased transcript abundance of carotenoid related genes including LRAT (LECITHIN RETINOL ACYLTRANSFERASE) and TS (TERPENE SYNTHASE) ( Figure 6 , Supplementary Data S1-9 ). β-carotene is cleaved into two molecules of retinal, which can then be reduced to retinol (vitamin A) (D’Ambrosio et al., 2011). LRAT plays a crucial role in metabolism of vitamin A by catalysing the formation of fatty acid retinyl esters from all-trans-retinol (Sears and Palczewski, 2016). The direct functional counterparts of LRAT are absent in plants, as plants do not require the synthesis of retinol. However, through a sequence-based homology analysis, homologues of LRAT have been identified in plants and are known as the H-Box/NC domain-containing gene family (Kaloudas and Penchovsky, 2018). The roles of LRAT homologs in plants still need investigation. However in Arabidopsis, NC domain containing proteins including AT5G06370, AT3G02700, and AT5G06370 regulate pollen germination and pollen tube growth (Wang et al., 2008). The higher abundance of LRAT homologs in Marchantia SL mutants aligns with their role in enhancing the germination of reproductive structures in Arabidopsis, which, in this study are gemmae. This suggests a conservation of certain functional domains across species, despite the differing metabolic requirements between plants and animals.

Plant terpene synthases (TS) are used to make monoterpenes, sesquiterpenes, and diterpenes that are then utilized in the production of ABA, chlorophyll and carotenoids (Cheng et al., 2007). TS are responsible for converting simple hydrocarbon skeletons into the vast array of terpenes important in plant defense, nodule formation, pollinator attraction, and interplant communication (Huang and Osbourn, 2019). For example, over expression of some terpenoid and terpene synthesis genes from Salvia officinalis in soybean are known to modulate rhizobia interaction and nodulation (Ali et al., 2021). The observed increase in LRAT and TS transcripts in SL precursor biosynthesis mutants may be attributed to either compensation for reduced metabolites controlled by SL precursor biosynthesis genes or because of a feedback mechanism where disrupted SL synthesis affects precursor synthesis. Additionally, Mpd27–2 had significantly higher levels of MpCCD7 ( Figure 7B ) suggesting that MpD27–2 and MpCCD7 both belong to the same metabolic pathway. When upstream components of a metabolic pathway are compromised then the amount of downstream enzyme could increase due to feedback regulation, compensation mechanisms, regulatory signals, or availability of substrate.

Alongside elevated expression of TS and LRAT, Mpd27–1, Mpd27–2 and Mpccd7 mutants exhibited an increased concentration of carotenoid-related metabolites ( Figure 8A ). Each of the three mutants showed a marked increase in the levels of cis-neoxanthin and lutein. The increased cis-neoxanthin could be attributed to downregulation of MpNCED in the SL mutants ( Supplementary Data S1-9 ). NCED plays a pivotal role in ABA biosynthesis, converting 9-cis-violaxanthin and 9′-cis-neoxanthin into xanthoxin, an ABA precursor (Perreau et al., 2020). Therefore, a decrease in MpNCED abundance could naturally result in the accumulation of neoxanthin, as less of it is converted into xanthoxin. Furthermore, the observed increase in lutein levels may also be attributed to disruptions in the activity of SL precursor biosynthesis enzymes, specifically those involved in cleaving 9-cis-zeaxanthin and 9-cis-lutein. For example, an in vitro study involving CCD7 enzymes sourced from rice, Arabidopsis, and garden pea demonstrated CCD7’s broad specificity. Notably, AtCCD7 and PsCCD7 were shown to effectively catalyze the conversion of 9-cis-lutein into 9-cis-3-OH-ϵ-apo-10′-carotenal (Baz, 2018). This evidence underscores CCD7’s role in carotenoid metabolism and its potential impact on lutein accumulation when the normal function of SL precursor biosynthesis enzymes is disrupted.

Mpd27–1 and Mpd27–2 mutants were observed to have elevated levels of β-carotene and α-carotene, respectively. This rise in carotenoid content within these mutants could be attributed to a decrease in isomerization processes, preventing their usual conversion into other compounds. The D27 enzyme from rice has been demonstrated in vitro to catalyze reactions across a variety of bicyclic carotenoid substrates, each containing at least one unsubstituted β-ionone ring, e.g., β-carotene, α-carotene, and cryptoxanthin (Bruno and Al-Babili, 2016). Similarly, AtD27 was capable of catalyzing the reverse isomerization of all-trans-/9-cis-β-carotene, showcasing enzyme activity that selectively targets specific isomeric forms (Abuauf et al., 2018). Notably, neither the rice nor the Arabidopsis D27 enzyme exhibits isomerization activity with 13-cis- or 15-cis-β-carotene, underscoring a distinct specificity for the C9-C10 double bond in these carotenoids. This specificity may play a critical role in the observed accumulation of β- and α-carotene in the d27 mutants, as it suggests a limited ability of the D27 enzyme to process these carotenoids into other forms, thereby contributing to their increased levels.

Mpd27–1 and Mpd27–2 both exhibited elevated levels of an unidentified carotenoid compound ( Figure 8A ). This observation suggests that MpD27–1 and MpD27–2 might play roles in the metabolism of carotenoid-like compounds, potentially beyond their established functions in carotene catabolism. Our hypothesis is further supported by a study where the enzymatic action of OsD27 on all-trans-α-carotene might be involved in producing a SL like compound such as heliolactone (Bruno and Al-Babili, 2016).

In this context, we propose that both MpD27–1 and MpD27–2, SL pathway isomerases have specific regulatory functions in the synthesis of SL compounds originating from β-carotene and α-carotene, respectively. This specific role of the two enzymes was evolutionarily favored to produce a diverse range of carotenoid-derived compounds, which might be implicated in controlling plant developmental processes such as regulation of release, germination, and growth of gemmae. Such a mechanism ensures that germination occurs under favorable growth conditions, which include sufficient ambient light and water availability.