Fay-Wei Li

Tomoaki Nishiyama

Manuel Waller

Eftychios Frangedakis

Jean Keller

Zheng Li

Noe Fernandez-Pozo

Michael S Barker

Tom Bennett

Miguel A Blázquez

Shifeng Cheng

Andrew C Cuming

Jan De Vries

Sophie De Vries

Pierre-Marc Delaux

Issa S Diop

Jill Harrison

Duncan Hauser

Jorge Hernández-García

Alexander Kirbis

John C Meeks

Isabel Monte

Sumanth K Mutte

Anna Neubauer

Dietmar Quandt

Tanner Robison

Masaki Shimamura

Stefan A Rensing

Juan Carlos Villarreal

Dolf Weijers

Susann Wicke

K.-S Gane

Wong

Keiko Sakakibara

Peter Szövényi
email: peter.szoevenyi@systbot.uzh.ch

Anthoceros genomes illuminate the origin of land plants and the unique biology of hornworts

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Land plants evolved from a Charophycean algal ancestor 470-515 million years ago [START_REF] Morris | The timescale of early land plant evolution[END_REF] and contributed to the greening of the terrestrial environment. The extant land plants consist of vascular plants and three bryophyte lineages-mosses, liverworts, and hornworts. While the phylogeny of land plants has been debated, recent evidence indicates that bryophytes are monophyletic with hornworts being sister to Setaphyta (liverworts plus mosses) [START_REF] De Sousa | Nuclear protein phylogenies support the monophyly of the three bryophyte groups (Bryophyta Schimp.)[END_REF][START_REF] Puttick | The Interrelationships of Land Plants and the Nature of the Ancestral Embryophyte[END_REF][START_REF] Wickett | Phylotranscriptomic analysis of the origin and early diversification of land plants[END_REF][START_REF]One Thousand Plant Transcriptomes Initiative. One thousand plant transcriptomes and phylogenomics of green plants[END_REF] .

The evolution of land plants is underlined by the rise of morphological, molecular and physiological innovations. Tracing the evolutionary origins of these key novelties is prone to errors due to the uncertainty in reconstructing the most recent common ancestor (MRCA) of land plants. More than 400 million years of independent evolution of the three bryophyte lineages have provided ample time for evolutionary changes to happen, and the limited availability of model systems for only two bryophyte lineages -mosses ( Physcomitrella patens ) [START_REF] Rensing | The Physcomitrella genome reveals evolutionary insights into the conquest of land by plants[END_REF] and liverworts ( Marchantia polymorpha ) [START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF] -makes the inferences even more difficult. Hornworts, as the earliest diverging lineage in bryophytes, are crucial to infer character evolution and reveal the nature of the MRCA of bryophytes and that of land plants.

Hornworts uniquely possess a combination of traits that connect them with both green algae and other land plant lineages [START_REF] Renzaglia | Comparative morphology and developmental anatomy of the Anthocerotophyta[END_REF] . For instance, most hornworts have a single chloroplast per cell with a pyrenoid capable of carrying out a carbon-concentrating mechanism (CCM) [START_REF] Smith | A pyrenoid-based carbon-concentrating mechanism is present in terrestrial bryophytes of the class Anthocerotae[END_REF] . Such pyrenoid-based CCMs cannot be found in any other land plants but frequently occur in algae [START_REF] Li | An overlooked window into carbon-concentrating mechanisms[END_REF] . Conversely, hornwort sporophytes are long-lived and moderately independent of gametophytes, which have been assumed to be a feature linking them to vascular plants [START_REF] Qiu | The deepest divergences in land plants inferred from phylogenomic evidence[END_REF] . Furthermore, hornwort sporophytes bear stomata which may be homologous with those of vascular plants [START_REF] Renzaglia | Hornwort stomata: architecture and fate shared with 400-Million-year-old fossil plants without leaves[END_REF] .

In addition to having characters exclusively shared with algae or with other land plants, hornworts also boast a wide range of distinctive biological features. For example, the presence of a basal sporophytic meristem and the asynchronous meiosis are unique to hornworts [START_REF] Renzaglia | New insights into morphology, anatomy, and systematics of hornworts[END_REF] . Moreover, hornworts are among the very few plants that have a symbiotic relationship with nitrogen-fixing cyanobacteria [START_REF] Adams | Cyanobacteria-bryophyte symbioses[END_REF] , and one particular hornwort species, Anthoceros punctatus , has been used as a model system to study plant-cyanobacteria interaction [START_REF] Meeks | Physiological adaptations in nitrogen-fixing Nostoc-plant symbiotic associations[END_REF] . Detailed genomic information on hornworts is essential to not only understand the evolutionary assembly of land plant-specific traits, but also to substantiate the full potential of hornworts as a model for studying the genetic basis of cyanobacterial symbiosis and pyrenoid-based CCM. Here we provide three high-quality genome assemblies and their annotations for the genus Anthoceros . We use these data to refine our inferences on the nature of land plant MRCA and to gain new insights into hornwort biology.

Results

Genome assembly and annotation

We assembled three hornwort genomes from Anthoceros agrestis (Bonn and Oxford strains) and A. punctatus . For A. agrestis Bonn, a combination of short-and long-read data with Chicago and Hi-C libraries resulted in a chromosomal-scale assembly with six largest scaffolds containing 95% of the assembled genome ( A. agrestis has 5-6 chromosome pairs; Fig. 1 and Supplementary Figure 1). For A. agrestis Oxford strain and A. punctatus , we used Oxford Nanopore sequencing to obtain high-quality assemblies composed of roughly 200 contigs with N50 over 1.7Mb (Table 1). The three genomes are highly collinear with a greater collinearity found between the two A. agrestis strains (Supplementary Figure 2; Supplementary Table 1). The collinearity, BUSCO scores (Supplementary Figure 3), and read mapping statistics (Supplementary Tables 23), show that the three genomes are of high quality and accuracy. The genome assemblies and annotations are available in CoGe ( https://genomevolution.org ; reviewer username: "hornwortgenome" and password: "lrm36ns89.DGcm"), and also on Figshare (private link: https://figshare.com/s/e3ebfc9104663c5d08de ).

The total assembly length varied between 117 and 133 Mb, which is consistent with the size estimates based on k -mer analysis (Table 1) but slightly larger than those from flow cytometry [START_REF] Szövényi | Establishment of Anthoceros agrestis as a model species for studying the biology of hornworts[END_REF][START_REF] Bainard | Genome size increases in recently diverged hornwort clades[END_REF] . Despite these genomes are among the smallest of land plants, their repetitive and transposable element contents are considerable (36-38%). Similar to other plant genomes, the most abundant repeats are LTR elements (>20%) followed by a large number of unclassified repeats and DNA elements. The genome size variation among the three strains can be largely attributed to the differences in repeat content (Supplementary Figure 4 and Supplementary Table 4). A combination of ab initio , evidence-based, and comparative gene prediction approaches resulted in 24,700-25,800 predicted protein-coding genes (Supplementary Table 5). The three hornwort genomes exhibit a high gene density compared to other land plants (Supplementary Table 6).

Anthoceros displays unusual centromere structure

The chromosomal-level assembly of A. agrestis Bonn revealed some peculiarities in the hornwort genome structures. In particular, we could not locate the typical vascular plant centromeric regions, which are usually composed of highly duplicated tandem repeats of 100-1000bp [START_REF] Jiang | A molecular view of plant centromeres[END_REF] . In A. agrestis Bonn, tandem repeats with a unit size over 30 bp gave rise to only very short arrays, and these repeats do not show a clear spatial clustering (Supplementary Figure 5). While gene density does fluctuate along the scaffolds, extensive regions with low gene density typical for centromeric regions of vascular plants were missing. Similarly, we could not identify stretches of scaffolds having an elevated repeat content (Fig. 1; Supplementary Figure 4), other than the putative telomeric regions. In other words, hornwort centromeres may not be characterized by a higher repeat density compared to other parts of the genome (see Supplementary Notes). Similar genome organizations were also discovered in the P. patens genome where genes and repeats are evenly distributed along the chromosomes [START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] . While it is tempting to hypothesize that this genomic organization may be a shared feature of bryophyte genomes, we nevertheless cannot rule out the possibility that the bona fide centromeres were not assembled properly despite the long-read and Hi-C data. Future work using immuno-labeling is necessary to confirm this hypothesis.

Limited collinearity across bryophyte and vascular plant genomes

A previous study on the moss P. patens genome implied that regions showing collinearity between the moss and some angiosperms may represent conserved collinear blocks since the MRCA of land plants [START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] . However, comparing bryophytes to vascular plants, shared ancestral gene blocks could not be identified, rather that the collinear regions with vascular plants were unique to each of the bryophyte genomes (Supplementary Figure 6; Supplementary Table 7). The largest number of genomic blocks collinear with at least one other land plant was found in the moss, followed by the liverwort and hornwort genomes (Supplementary Figure 6). Within bryophytes, no collinear segment conserved across all the three lineages was found, although there were genomic regions exclusively collinear between each of two bryophyte genomes (Supplementary Figure 6). In general, there was more collinearity between the liverwort and the moss than between the hornwort and the liverwort/moss genomes. The number of such collinear regions, however, were small compared to those detected across vascular plants (Supplementary Figure 6). Altogether, these findings imply that the deep divergence of the moss, hornwort, and liverwort genomes may have led to limited collinearity both among bryophytes, as well as between bryophytes and vascular plants. Our results also suggest that each bryophyte genome potentially retained a unique, non-overlapping set of collinear regions from the MRCA of land plants.

Absence of large scale genome duplication in Anthoceros

Whole genome duplications (WGD) have played an important role in shaping plant evolution and possibly underlie several adaptive radiations [START_REF] Landis | Impact of whole-genome duplication events on diversification rates in angiosperms[END_REF] . A previous study, based on Ks divergence in transcriptomic datasets, suggested that hornworts may not have experienced any WGD event [START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] , similar to M. polymorpha [START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF] and Selaginella moellendorffii [START_REF] Banks | The Selaginella genome identifies genetic changes associated with the evolution of vascular plants[END_REF] . Our Ks plots on the annotated Anthoceros genes similarly show no sign of WGD (Supplementary Figure 7). To further corroborate this, we investigated patterns of intra-genomic synteny in the three hornwort genomes, as well as the published M. polymorpha and P. patens genomes for comparison. We found very little self-synteny in the hornwort genomes (Fig. 1b-d), thus providing a strong evidence for the lack of WGD in Anthoceros . The high proportion of the genomes that are not syntenic is comparable to that in M. polymorpha (Fig. 1e). On the other hand, P. patens shows a clear 1:1 (and some 1:2) self-syntenic relationship (Fig. 1f), which is consistent with the earlier report and indicative of two rounds of WGD [START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] .

Small repertoire of transcription associated proteins

We found that 2.4-2.6% of the gene set of the three Anthoceros genomes was annotated as transcription associated proteins (TAPs) (Fig. 2; Supplementary Table 8). Compared to other land plants [START_REF] Wilhelmsson | Comprehensive genome-wide classification reveals that many plant-specific transcription factors evolved in streptophyte algae[END_REF] , this is on the very low end of the spectrum. Furthermore, approximately 2/3rd (56) of the hornwort TAP families were smaller in size than in M. polymorpha . Given such a minimal TAP repertoire, hornworts can serve as an excellent baseline model to study the evolution and diversification of transcriptional networks. Despite its streamlined nature, some TAPs were only found in hornworts and vascular plants but not in the other two bryophyte genomes, with the most notable example being YABBY (Supplementary Table 8; Supplementary Note). Such TAPs likely evolved in the MRCA of land plants but were lost in the mosses and liverworts. We also detected TAP families that were present in all streptophytes but lost either in the hornwort genomes (e.g. SRS transcription factor [TF]) or in M. polymorpha (e.g. type I MADS-box TF). Altogether, our findings suggest a dynamic TAP family turnover in the early evolution of land plants with multiple independent losses in different bryophyte lineages.

Genes related to sporophyte development

While hornworts have a gametophyte-dominant life cycle like other bryophytes, their sporophyte generation exhibits several unique features [START_REF] Ligrone | The origin of the sporophyte shoot in land plants: a bryological perspective[END_REF] . First, after fertilization the zygote division in hornworts is longitudinal, whereas zygotes in all other land plants undergo transverse division. Second, the hornwort sporophyte maintains a basal sporophytic meristem that continuously differentiates into mature tissues towards the tip. A common origin of the indeterminate sporophyte development in hornworts and the vascular plant shoot apical meristem (SAM) has been hypothesized [START_REF] Ligrone | The origin of the sporophyte shoot in land plants: a bryological perspective[END_REF] . Lastly, hornwort sporophytes have stomata (Fig. 2c) similar to mosses and vascular plants, and the basic regulation may be shared across all stomatous lineages of land plants [START_REF] Chater | Origins and Evolution of Stomatal Development[END_REF] . Nevertheless, firm evidence supporting the homology of meristems as well as stomata is scarce. Here we found that multiple genes critical for flowering plant SAM and stomata function have homologs in the hornwort genomes and are preferentially expressed in the sporophyte phase.

Class 1 Knotted1 -like homeobox (KNOX1) genes regulate sporophytic meristem activity in both P. patens and vascular plants [START_REF] Hay | KNOX genes: versatile regulators of plant development and diversity[END_REF] , while Class 2 Knotted1 -like homeobox (KNOX2) genes maintain the sporophyte cell fate in P. patens [START_REF] Sakakibara | KNOX2 genes regulate the haploid-to-diploid morphological transition in land plants[END_REF] . Interestingly, the KNOX1 ortholog is lost in the Anthoceros genomes, and only KNOX2 genes were found (Supplementary Figure 8; Supplementary Tables 89). The KNOX2 orthologs showed a strong sporophyte-specific expression (Fig. 2), which implies that KNOX2's involvement in maintaining sporophytic cell fate may be conserved in all land plants. Heterodimerization of KNOX1/KNOX2 and BELL-LIKE HOMEOBOX proteins is a deeply conserved molecular mechanism that is required for the KNOX functions [START_REF] Arun | Convergent recruitment of TALE homeodomain life cycle regulators to direct sporophyte development in land plants and brown algae[END_REF] . We found that in hornworts, a single BELL and a single KNOX2 gene were expressed in the early stages of sporophyte development (Fig. 2; Supplementary Tables 89). This suggests that hornwort sporophyte identity may be determined by KNOX2 through interaction with BELL .

WUSCHEL-related homeobox 13 like (WOX13L) genes are involved in the zygote development and stem cell formation in the moss P. patens [START_REF] Sakakibara | WOX13-like genes are required for reprogramming of leaf and protoplast cells into stem cells in the moss Physcomitrella patens[END_REF] . A. thaliana WOX13 promotes replum formation in the fruit [START_REF] Romera-Branchat | The WOX 13 homeobox gene promotes replum formation in the Arabidopsis thaliana fruit[END_REF] and WOX14 promotes vascular cell differentiation [START_REF] Denis | WOX14 promotes bioactive gibberellin synthesis and vascular cell differentiation in Arabidopsis[END_REF] .

The Anthoceros genomes have four WOX13L members (Supplementary Figure 8; Supplementary Tables 89), and WOX13La is specifically expressed in sporophytes while WOX13Lbcd have expression at both gametophyte and sporophyte generations (Fig. 2b) and may have diverse roles in stem cell maintenance and sporophyte development. The Anthoceros genomes also have a single FLORICAULA / LEAFY ( FLO / LFY ) gene (Supplementary Figure 8; Supplementary Tables 89), which in P. patens and A. thaliana controls zygote development and SAM maintenance, respectively [START_REF] Tanahashi | Diversification of gene function: homologs of the floral regulator FLO/LFY control the first zygotic cell division in the moss Physcomitrella patens[END_REF] . In hornworts, LFY is predominantly expressed in the gametophyte stages (Fig. 2) while in P. patens it is expressed both in the gametophyte and the sporophyte. It is possible that such differences may contribute to the unique developmental pattern of hornwort sporophytes.

Stomatal development in A. thaliana and P. patens is regulated by a conserved genetic toolbox, including the basic helix-loop-helix (bHLH) transcription factors SMF ( SPCH , MUTE , and FAMA ), ICE/SCREAMs ( SCRM s), EPIDERMAL PATTERNING FACTOR ( EPF ), ERECTA , and TOO MANY MOUTHS ( TMM ) genes [START_REF] Lee | The plant stomatal lineage at a glance[END_REF][START_REF] Chater | Origin and function of stomata in the moss Physcomitrella patens[END_REF] . FAMA in particular is involved in the final guard cell differentiation and serves as the key switch. Orthologs of SMF , TMM, and EPF were absent in M. polymorpha , consistent with the fact that liverworts do not have stomata [START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF] . We found orthologs of FAMA (SMF) , SCRM , ERECTA , EPF , and TMM in the Anthoceros genomes (in line with a previous study based on our earlier genome draft [START_REF] Chater | Origins and Evolution of Stomatal Development[END_REF] ; Supplementary Table 10, Supplementary Figure 9). SMF , SCRM , TMM and EPF showed sporophyte-specific expression pattern (Fig. 2), suggesting that they may have similar roles in stomata patterning in hornworts. While ERECTA was also expressed during early sporophyte development, its expression fluctuated between replicates and were inconclusive. EPF expression showed similar inconsistency among replicates, but did not influence our conclusion about its sporophyte specific expression. In addition to EPF , an EPF -like gene in the EPFL4-6 clade was found in hornworts (Supplementary Figure 9), which is specifically expressed in gametophytes with a higher expression toward maturity and thus perhaps involved in a different cell-cell signaling other than stomatal regulation. EPF4 and EPF6 in A. thaliana are involved in coordination of central and peripheral zone in shoot apical meristem [START_REF] Kosentka | EPFL Signals in the Boundary Region of the SAM Restrict Its Size and Promote Leaf Initiation[END_REF] . Taken together, our data are consistent with a single origin of stomatal differentiation mechanism among all stomatous land plants, though the positional determination may have evolved differently (Supplementary Notes).

Genes related to phytohormone synthesis and signaling

The Anthoceros genomes contain the genetic chassis for the biosynthesis and signaling of abscisic acid, auxin, cytokinin, ethylene, and jasmonate (see Supplementary Notes; Supplementary Figures 101112; Supplementary Table 10), reaffirming the origins of these pathways in the MRCA of land plants [START_REF] Rensing | The Physcomitrella genome reveals evolutionary insights into the conquest of land by plants[END_REF][START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF][START_REF] Nishiyama | The Chara genome: secondary complexity and implications for plant terrestrialization[END_REF] . Similar to M. polymorpha and P. patens , salicylic acid signaling components are found in hornworts, but not the receptor-related genes. While DELLA is present, orthologs of gibberellin (GA) receptor GID1 and GA oxidases are missing from the Anthoceros genomes. This is consistent with the recent hypothesis that DELLA was recruited to GA signaling pathway later in plant evolution [START_REF] Hernandez-Garcia | Origin of gibberellin-dependent transcriptional regulation by molecular exploitation of a transactivation domain in DELLA proteins[END_REF] . Hornworts also possess enzymes to synthesize strigolactones but genes involved in strigolactone signaling are absent. This supports the idea that strigolactones are an ancient non-hormonal signal for rhizospheric communication with mycorrhizal fungi 37

Genetic network for arbuscular mycorrhizal symbiosis was present in the MRCA of land plants

The symbiotic relationship with arbuscular mycorrhizal fungi (AMF) is one of the key innovations underlying plants' successful colonization and diversification on land. Evidence of AMF can be found inside plant megafossils 407 million years ago [START_REF] Remy | Four hundred-million-year-old vesicular arbuscular mycorrhizae[END_REF][START_REF] Strullu-Derrien | Fossil filamentous microorganisms associated with plants in early terrestrial environments[END_REF] , and in almost all extant plant lineages (hornworts, liverworts, and vascular plants). Recent genetic studies have identified a suite of genes in the angiosperms that regulate the establishment and maintenance of AMF symbiosis. Some of these genes are also required for legume-rhizobial interaction, and are often referred to as the common symbiosis genes [START_REF] Parniske | Arbuscular mycorrhiza: the mother of plant root endosymbioses[END_REF] .

While a few components can be traced back to as far as charophyte algae [START_REF] Delaux | Algal ancestor of land plants was preadapted for symbiosis[END_REF] , the question of when exactly did the entire AM symbiosis genetic network originate still remains open. This is partly because both of the bryophytes that have published genomes to date (i.e. P. patens and M. polymorpha ) are incapable of AMF symbiosis and may have secondarily lost the symbiosis genes, as exemplified in some angiosperms [START_REF] Delaux | Comparative phylogenomics uncovers the impact of symbiotic associations on host genome evolution[END_REF] . Here we show that all the key angiosperm AMF symbiosis genes have orthologs in the three hornwort genomes (Fig. 3, Supplementary Table 11, Supplementary Figure 13). Although their roles in hornwort AM symbiosis remains to be tested, this result provides strong evidence that the genetic infrastructure required for AM symbiosis was already present in the MRCA of land plants. Importantly, the presence of these genes in liverworts [START_REF] Delaux | Algal ancestor of land plants was preadapted for symbiosis[END_REF] and hornworts makes this conclusion insensitive to any uncertainty of the land plant phylogeny. We have not succeeded in reconstituting hornwort-AMF symbiosis in vitro and hence are unable to test these ortholog expressions in the context of AMF. Nevertheless, we found that in both A. agrestis (Oxford strain) and A. punctatus , one of the AMF symbiosis genes, RAM1 , was up-regulated when plants were nitrogen-starved (Fig. 3). Nitrogen limitation is a major trigger for cyanobacteria symbiosis in hornworts, which might implicate RAM1's involvement in symbiosis but further genetic studies are needed.

Genes related to cyanobacterial symbiosis

Symbiosis with nitrogen-fixing cyanobacteria is a rare trait, with limited appearances in a few plant lineages [START_REF] Adams | Cyanobacteria-bryophyte symbioses[END_REF] . In bryophytes, although mosses frequently harbor epiphytic cyanobacteria [START_REF] Rousk | Moss-cyanobacteria associations as biogenic sources of nitrogen in boreal forest ecosystems[END_REF] , only hornworts and two liverwort species host cyanobacteria endophytically within specialized slime-filled cavities [START_REF] Adams | Cyanobacteria-bryophyte symbioses[END_REF] . Amongst all the plant associations with cyanobacteria, most of the research has been done on hornworts, using A. punctatus (sequenced here) and the cyanobacterium Nostoc punctiforme as the study system.

Although several cyanobacterial genes from N. punctiforme have been identified that are key to initiation of symbiotic association [START_REF] Meeks | Physiological adaptations in nitrogen-fixing Nostoc-plant symbiotic associations[END_REF] , nothing is known about the hornwort genetics. Here we generated RNA-seq data to compare the gene expression of symbiont-free (either nitrogen-starved or nitrogen-fed) and symbiosis-reconstituted hornworts (Fig. 3). This experiment was conducted with both A. punctatus and the A. agrestis Oxford isolate. We identified 40 genes that, when the cyanobionts are present, are highly induced (>16-fold) in both hornwort species (Fig. 3; Supplementary Table 12). These include a number of receptor kinases, transcription factors, and transporters. Of particular interest is a SWEET sugar transporter in the SWEET16/17 clade (Fig. 3; Supplementary Figure 14), which is minimally transcribed under the symbiont-free states but is among the highest expressed genes in symbiosis (>10 3 fold-change). The up-regulation of SWEET in symbiosis is interesting because it implies that this sugar transporter is dedicated to supplying carbon rewards to the cyanobionts. This implication is supported by the fact that only exogenous glucose, fructose or sucrose sustained dark nitrogen fixation in the A . punctatus -N . punctiforme association [START_REF] Steinberg | Physiological sources of reductant for nitrogen-fixation activity in Nostoc sp. strain UCD 7801 in symbiotic association with Anthoceros punctatus[END_REF] , and the observation that inactivation of a carbohydrate permease in N . punctiforme resulted in a defective symbiotic phenotype [START_REF] Ekman | A Nostoc punctiforme sugar transporter necessary to establish a cyanobacterium-plant symbiosis[END_REF] . Parallelly, SWEET is involved in mycorrhizal symbiosis as well, but a different ortholog, in the SWEET1 clade, was recruited [START_REF] An | A Medicago truncatula SWEET transporter implicated in arbuscule maintenance during arbuscular mycorrhizal symbiosis[END_REF] .

Another gene of interest is subtilase. Members of this gene family have been shown to be highly up-regulated in a wide variety of microbial symbioses, including rhizobial [START_REF] Kistner | Seven Lotus japonicus genes required for transcriptional reprogramming of the root during fungal and bacterial symbiosis[END_REF] , mycorrhizal [START_REF] Takeda | Apoplastic plant subtilases support arbuscular mycorrhiza development in Lotus japonicus[END_REF] , and actinorhizal [START_REF] Fournier | Cell remodeling and subtilase gene expression in the actinorhizal plant Discaria trinervis highlight host orchestration of intercellular Frankia colonization[END_REF][START_REF] Ribeiro | A nodule-specific gene encoding a subtilisin-like protease is expressed in early stages of actinorhizal nodule development[END_REF][START_REF] Svistoonoff | cg12 expression is specifically linked to infection of root hairs and cortical cells during Casuarina glauca and Allocasuarina verticillata actinorhizal nodule development[END_REF] interactions. RNAi knockdown of a subtilase ( SBTM1 ) in the legume Lotus japonicus also resulted in a decreased arbuscule formation [START_REF] Takeda | Apoplastic plant subtilases support arbuscular mycorrhiza development in Lotus japonicus[END_REF] . Here we found that in both A. punctatus and A. agrestis , a subtilase homolog was similarly induced by cyanobacteria symbiosis. Phylogenetic reconstruction showed that this hornwort subtilase is not orthologous to those involved in other plant symbioses (Supplementary Figure 15). Taken together, our results imply that hornworts might have convergently recruited SWEET and subtilase for cyanobacterial symbiosis, although in both cases not the same orthologs were used as in other plant-microbe symbioses.

Pyrenoid-based carbon concentrating mechanism

To enable a more efficient photosynthesis, hornworts, cyanobacteria, and many eukaryotic algae have evolved biophysical carbon-concentrating mechanism (CCM) inside individual chloroplasts [START_REF] Meyer | The algal pyrenoid: key unanswered questions[END_REF] . In this case, chloroplasts use inorganic carbon transporters and carbonic anhydrases to locally concentrate CO 2 in the pyrenoids, a specialized chloroplast compartment where RuBisCOs aggregates. Pyrenoids can thus boost photosynthetic efficiency and reduce photorespiration. Such pyrenoid-based CCM has been extensively studied in the model green alga Chlamydomonas reinhardtii with the hope to install a CCM in crop plants [START_REF] Rae | Progress and challenges of engineering a biophysical CO2-concentrating mechanism into higher plants[END_REF] .

Hornworts are the only land plants with a pyrenoid-based CCM. Interestingly, for the past 100 million years, pyrenoids in hornworts are inferred to have been repeatedly lost and gained [START_REF] Villarreal Aguilar | Hornwort pyrenoids, carbon-concentrating structures, evolved and were lost at least five times during the last 100 million years[END_REF] , suggesting that pyrenoid development and function is controlled by a few master switches. The genetics behind hornwort pyrenoids, however, has remained completely unknown. We hence explored whether hornwort genomes have genes that are known to be required for pyrenoid-based CCM in C. reinhardtii . While many of the C. reinhardtii CCM genes [START_REF] Rae | Progress and challenges of engineering a biophysical CO2-concentrating mechanism into higher plants[END_REF] do not have clear homologs in hornworts (nor in any other land plants), we did find LCIB (low-CO 2 inducible B) to be present in the hornwort genomes and 1KP transcriptomes [START_REF]One Thousand Plant Transcriptomes Initiative. One thousand plant transcriptomes and phylogenomics of green plants[END_REF] (Fig. 4). Apart from hornworts, no LCIB homolog could be found in other plant genomes sequenced to date. The uniquely shared presence of LCIB in pyrenoid-bearing algae and hornworts implies that LCIB might have a role in the hornwort CCM. The phylogenetic tree indicates that the hornwort LCIBs form a clade sister to the Klebsormidium nitens homolog (Fig. 4) and thus is consistent with the organisms tree with many losses in various lineages. In this scenario, the MRCA of land plants had LCIB.

In C. reinhardtii , LCIB gene expression is highly induced by CO 2 limitation, and the encoded proteins localize around pyrenoids to presumably block CO 2 leakage [START_REF] Wang | LCIB in the Chlamydomonas CO2-concentrating mechanism[END_REF][START_REF] Atkinson | Introducing an algal carbon-concentrating mechanism into higher plants: location and incorporation of key components[END_REF] . All the hornwort LCIB sequences have the conserved amino acid residues at the active sites that are shared with other algal LCIBs [START_REF] Jin | Structural insights into the LCIB protein family reveals a new group of β-carbonic anhydrases[END_REF] (Fig. 4b). However, unlike C. reinhardtii , we did not find LCIB to be differentially expressed when plants are grown at different CO 2 levels (Supplementary Figure 12). This, nevertheless, cannot rule out LCIB 's involvement in CCM because hornwort CCM was reported to be constitutively expressed and not regulated by CO 2 level [START_REF] Hanson | Diffusion Limitation and CO2 Concentrating Mechanisms in Bryophytes[END_REF] . Whether LCIB homologs have a similar function and localization in hornworts remains to be experimentally tested.

Discussion

The hornwort genomes presented here offer a unique window into the biology of land plant MRCA. For example, the Anthoceros genomes lack KNOX1 , while P. patens and M. polymorpha lack YABBY genes. This suggests that the MRCA of land plants had both of these key developmental genes, and independent gene losses occurred in different bryophyte lineages. While LEAFY expression is predominantly in the gametophyte stage, YABBY , KNOX2 , BELL and some WOX13L genes are up-regulated in the hornwort sporophytes (Fig. 2). In addition, several stomata related genes are present in the Anthoceros genomes and expressed in early sporophyte development (Fig. 2), implying a homology of stomata at the genetic level. Finally, we found that the genes required for AM symbiosis are conserved in Anthoceros (Fig. 3), providing the first solid evidence that the MRCA of land plants was already equipped with the genetic network for AM symbiosis. In-depth analysis on the evolution of the plant hormones (abscisic acid, auxin, gibberellin, jasmonate, salicylic acid, and strigolactone), light signaling, peptidoglycan synthesis, and chloroplast development can be found in Supplementary Notes.

The Anthoceros genomes shared several features with the two other published bryophyte genomes. Most notably is the absence of tandem repeats that make up the typical centromeric regions. Further studies are needed to identify the centromeric regions and understand its structure. While P. patens has experienced two rounds of WGD [START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] , none can be found in Anthoceros and M. polymorpha (Fig. 1). This might explain the minimal representation of transcription factors in the latter two genomes. Furthermore, our functional genomic data shed light on the genetic framework that underpins features that are unique to hornworts. We identified a suite of candidate genes underlying hornwort-cyanobacteria symbiosis (Fig. 3). This includes a SWEET transporter that might be involved in nutrient transfer with the cyanobionts. A well-characterized C. reinhardtii CCM gene, LCIB , was conserved in hornworts but apparently lost in all other plant lineages (Fig. 4). Whether or not LCIB also participates in hornwort CCM awaits future functional characterization.

The recent advances of "seed-free genomics" have significantly improved our understanding of streptophyte evolution [START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF][START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF][START_REF] Banks | The Selaginella genome identifies genetic changes associated with the evolution of vascular plants[END_REF][START_REF] Nishiyama | The Chara genome: secondary complexity and implications for plant terrestrialization[END_REF][START_REF] Li | Fern genomes elucidate land plant evolution and cyanobacterial symbioses[END_REF][START_REF] Hori | Klebsormidium flaccidum genome reveals primary factors for plant terrestrial adaptation[END_REF] . Here our hornwort genomes fill in yet another critical gap, and are beginning to illuminate the dawn of land plants as well as the unique biology of hornworts.

Methods

Plant materials

Cultures of Anthoceros agrestis

(Oxford and Bonn strains) and A. punctatus were all derived from a single spore, and axenically propagated and maintained on either BCD 61 or Hatcher's 62 medium. Supplementary Table 13 shows the origin and specimen voucher for each of the three strains.

Chromosome count

The tip of an A. agrestis Oxford thallus was cut into small pieces and fixed with 4% glutaraldehyde in 0.05M phosphate buffer (pH 7.0) for 12 hours at 4°C. After washing with the buffer for 10 minutes, cell walls were digested for 2 hours with a solution containing 1% Driselase (Sigma-Aldrich, St. Louis, MO, USA), 1% Cellulase Onozuka RS (Yalult, Tokyo, Japan), 1% Pectolyase (Kikkoman, Tokyo, Japan), 0.5% IGEPAL CA-630, and 1% bovine serum albumin (BSA) at 30°C. After several washes with the buffer, the samples were incubated in 0.05M phosphate buffer containing 0.1% TritonX-100 for 12 hours at 4°C. After several washes with the buffer, the samples were transferred onto MAS coated slide glasses (Matsunami Glass, Osaka, Japan) and coverslipped. The slides were then pressed with a thumb directly over the coverslip. After removal of the coverslip, the slides were air-dried for 10 minutes at room temperature and then extracted with methanol at -20°C for 10 min. After the staining with the buffer containing 1 μg/L 4,6-diamidino-2-phenylindole (DAPI) for 5 minutes, the slides were mounted with Vectashield mounting medium (Vector Laboratory Burlingame, CA, USA) and observed with a fluorescence microscope under UV-light excitation.

DNA sequencing

Hornwort DNA was extracted using a CTAB-precipitation method modified from [START_REF] Nagar | High purity, high molecular weight DNA extraction from rust spores via CTAB based DNA precipitation for long read sequencing v1[END_REF] . Nanopore libraries were prepared by SQK-LSK108, and sequenced on MinION R9 flow cells for 48 hours. Basecalling was done by Albacore.

For A. agrestis Bonn, the TrueSeq DNA Nano kit (Illumina) was used to prepare sequencing libraries which were sequenced (PE 150 bp) on HiSeq4000 at the Functional Genomic Center Zurich (FGCZ). For A. agrestis Oxford 251 PE reads, a PCR-Free library was prepared using a KAPA Hyper prep kit according to the protocol published by Broad Institute 64 ). The library was mixed (5%) with other barcoded libraries, and sequenced on Illumina HiSeq1500 (2 lanes with Rapid mode; OnBoardClustering) at the National Institute of Basic Biology (NIBB). For A. punctatus , Illumina genomic libraries were prepared by BGI and sequenced on HiSeq4000. Read quality and adaptor trimming was done by fastp [START_REF] Chen | fastp: an ultra-fast all-in-one FASTQ preprocessor[END_REF] with the default setting. For A. agrestis Bonn, additional Chicago and HiC libraries were prepared by DoveTail Genomics. A total of two Chicago libraries and one HiC library were prepared with a physical coverage of 300 and 200x.

To calculate the read mapping rates, trimmed reads were mapped to the final assemblies using bwa mem [START_REF] Li | Fast and accurate short read alignment with Burrows-Wheeler transform[END_REF] . The mean and median insert sizes were calculated using picard collectInsertSizemetrics. Unmapped reads were counted with samtools view -f 4 [START_REF] Li | The sequence alignment/map format and SAMtools[END_REF] and divided with the total number of reads to calculate percent mapped. Reads mapped to chloroplast (Cp) and mitochondrial (Mt) genomes were counted with samtools and divided by the number of records in the bam file that is not flagged as unmapped to obtain percentage of Cp and Mt reads. All the raw sequences are deposited to NCBI SRA under the BioProject PRJNA574424, PRJNA574453, and to ENA under the study accessions PRJEB34763, PRJEB34743 (Supplementary Tables 23).

Genome assembly

Genome sizes for the three Anthoceros were first estimated by Jellyfish [START_REF] Marçais | A fast, lock-free approach for efficient parallel counting of occurrences of k-mers[END_REF] in conjunction with GenomeScope [START_REF] Vurture | GenomeScope: fast reference-free genome profiling from short reads[END_REF] . Draft assembly for A. agrestis Bonn strain was first generated using a hybrid approach including Oxford nanopore (approx. 60x) and Illumina paired-end reads (approx. 150x) using MaSuRCA version 3.2.8 [START_REF] Zimin | The MaSuRCA genome assembler[END_REF] . After assembly, base call quality was improved by two rounds of pilon polishing [START_REF] Walker | An integrated tool for comprehensive microbial variant detection and genome assembly improvement[END_REF] . We mapped Chicago and Hi-C reads back to the draft assembly and used DoveTail's HiRise assembler v2.1.2 [START_REF] Putnam | Chromosome-scale shotgun assembly using an in vitro method for long-range linkage[END_REF] for scaffolding. Contigs of the draft assembly were first scaffolded using the Chicago library to correct smaller scale errors and improve contiguity. Finally, the output assembly was further scaffolded using the HiC libraries and DoveTail's HiRise assembler v2.1.2 [START_REF] Putnam | Chromosome-scale shotgun assembly using an in vitro method for long-range linkage[END_REF] to derive the final assembly.

Genome assemblies of A. agrestis Oxford strain and A. punctatus were generated with the minimap2-miniasm assembler [START_REF] Li | Minimap and miniasm: fast mapping and de novo assembly for noisy long sequences[END_REF] using only the nanopore reads. We then used four iterations of minimap2-racon [START_REF] Vaser | Fast and accurate de novo genome assembly from long uncorrected reads[END_REF] to derive the consensus sequence, followed by six rounds pilon polishing [START_REF] Walker | An integrated tool for comprehensive microbial variant detection and genome assembly improvement[END_REF] .

Contamination removal

While our cultures were grown in a putative axenic condition, low level of contamination cannot be completely ruled out. We therefore used blobtools [START_REF] Laetsch | Interrogation of genome assemblies[END_REF] to identify scaffolds/contigs primarily consisting of contaminant sequences. The Hi-C library theoretically should sort DNA sequences originating from different organisms because cross-linking occurred within the nuclei. Therefore, we hypothesized that dropping scaffolds mainly with non-Streptophyte affiliation will effectively remove contaminants from our assembly. For A. agrestis Bonn, we used both the full uniprot and the NCBI nt database and blobtools to assign the taxonomic affiliation to each scaffold with an e-value of 10 -4 . We found that some of the small scaffolds were classified as Ascomycota and Cyanobacteria; these scaffolds were then removed from the assembly. For the A. agrestis Oxford and A. punctatus genomes, assemblies were contamination-filtered in a similar way.

RNA-seq dataset and analysis on developmental stages

To study the expression pattern of transcription factor genes across developmental stages, we generated RNA-seq libraries for the following stages of the A. agrestis Bonn strain in two biological replicates: (1) spores after two weeks of germination, (2) four-week-old gametophytes, (3) two-month-old gametophytes, (4) sporophytes shorter than 5 mm, (5) sporophytes 5-10 mm, (6) sporophytes longer than 1 cm with brown or black tips. Plants were grown on agar plates containing BCD medium 61 at 22°C. RNA was extracted with the Spectrum Total RNA Plant Kit (Sigma-Aldrich) and stranded RNA-seq libraries were prepared using the TrueSeq Stranded mRNA Library Prep kit (Illumina). Libraries were sequenced at the FGCZ on a HiSeq4000 machine. We used trimmomatic [START_REF] Bolger | Trimmomatic: a flexible trimmer for Illumina sequence data[END_REF] to quality filter and trim the raw reads. Gene expression was estimated using Salmon [START_REF] Patro | Salmon provides fast and bias-aware quantification of transcript expression[END_REF] and differential expression done by DESeq2 (log2fold ≥ 2, false-discovery rate ≤ 0.05, and normalized reads counts) [START_REF] Love | Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2[END_REF] .

We also generated separate thallus RNA-seq data for the Oxford strain (for annotation purpose). The plants were cultured on solid BCD plates, and total RNA was extracted using RNeasy Plant Mini kit (QIAGEN). The library was prepared using the TrueSeq stranded mRNA Library Prep kit (Illumina) and sequenced on Hiseq1500.

RNA-seq dataset and analysis on cyanobacterial symbiosis

Liquid cultures of A. agrestis Oxford and

A. punctatus were used in this experiment. To establish liquid cultures, plants were transferred from solid BCD plates to flasks with 100 ml BCD media solution, and placed on an orbital shaker with 130 rpm for two weeks. For the cyano-/N+ and cyano-/N-conditions, plants were transferred to fresh new BCD solution with and without KNO 3 , respectively and grown for 10 days before harvest. To reconstitute cyanobacterial symbiosis (with Nostoc punctiforme ATCC 29133), we followed the method of Enderlin and Meeks [START_REF] Enderlin | Pure culture and reconstitution of the Anthoceros-Nostoc symbiotic association[END_REF] , but using BCD as the growth medium. Three biological replicates were done for each condition. RNA was extracted by the Spectrum Total RNA Plant Kit (Sigma-Aldrich). The Illumina libraries were prepared by BGI and sequenced on HiSeq4000. Sequencing reads were mapped to the respective genomes by HiSat2 [START_REF] Kim | HISAT: a fast spliced aligner with low memory requirements[END_REF] , and transcript abundance quantified by Stringtie [START_REF] Pertea | StringTie enables improved reconstruction of a transcriptome from RNA-seq reads[END_REF] . We used DESeq2 [START_REF] Love | Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2[END_REF] to carry out differential gene expression analysis, with false discovery rate set to 0.005 and log2fold change threshold set to 1. To identify genes that are differentially expressed in both A. agrestis Oxford and A. punctatus , we used on the Orthofinder gene family classification results (see below) coupled with phylogenetic analysis if needed.

RNA-seq dataset and analysis on CO 2 response

For the CO 2 experiment, we grew hornworts in magenta boxes with vented lids in order to allow air circulation while maintaining sterility. A. agrestis Oxford strain was used in this experiment and kept on solid BCD medium. We subjected the plant cultures to one of the three CO 2 environments at 150 (low), 400 (ambient), and 800 (high) ppm in a CO 2 -controlled growth chamber for 10 days (12/12hr day/night cycle). Three biological replicates were done for each treatment. RNA was extracted by the Spectrum Total RNA Plant Kit (Sigma-Aldrich). The Illumina libraries were prepared by BGI, and sequenced on HiSeq4000. One of the low CO 2 samples failed to produce high quality library, and as a result the low CO 2 condition has only two replicates. RNA-seq data analysis was done following the same procedure as described above. We used BiNGO [START_REF] Maere | a Cytoscape plugin to assess overrepresentation of gene ontology categories in biological networks[END_REF] for gene ontology enrichment analysis, and REVIGO [START_REF] Supek | REVIGO summarizes and visualizes long lists of gene ontology terms[END_REF] to summarize and visualize the results.

Repeat annotation

For repeat annotation, we first built custom repeat libraries for each genome using RepeatModeler [START_REF] Smit | RepeatModeler Open-1.0[END_REF] and LTR_retriever [START_REF] Ou | LTR_retriever: A highly accurate and sensitive program for identification of long terminal repeat retrotransposons[END_REF] . The libraries were filtered to remove protein-coding genes by blasting against the UniProt plant database. We then used RepeatMasker [START_REF] Smit | RepeatMasker Open-4.0[END_REF] to annotate and mask the repetitive regions for each genome.

RNA-seq, transcript, and protein evidence

We pooled A. agrestis Bonn, Oxford and A. punctatus RNA-seq reads together and mapped them onto each of the genome assemblies using HiSat2 [START_REF] Kim | HISAT: a fast spliced aligner with low memory requirements[END_REF] . We used all RNA-seq evidence available owing to the low nucleotide divergence among the three genomes. Transcriptomes were assembled for each species/strain separately. We used Portcullis 87 to filter out bad splice junctions, and Stringtie 81 to assemble the transcripts. We additionally used Trinity [START_REF] Grabherr | Full-length transcriptome assembly from RNA-Seq data without a reference genome[END_REF] to generate both de novo and genome-guided transcriptome assemblies. We combined Trinity transcripts using the PASA pipeline [START_REF] Haas | Improving the Arabidopsis genome annotation using maximal transcript alignment assemblies[END_REF] and derived high-quality transcripts with Mikado [START_REF] Venturini | Leveraging multiple transcriptome assembly methods for improved gene structure annotation[END_REF] . To obtain protein homology information, we retrieved the 19 proteomes (only primary transcripts; Supplementary Table 14) and aligned them to the genome assemblies using exonerate [START_REF] Slater | Automated generation of heuristics for biological sequence comparison[END_REF] . We kept only hits with at least 60% coverage and a similarity above 60%.

Gene prediction

We used RNA-seq, transcript, and protein evidence to train Augustus 92 within Braker2 [START_REF] Hoff | BRAKER1: Unsupervised RNA-seq-based genome annotation with GeneMark-ET and AUGUSTUS[END_REF] . Because the resulting gene models were heavily dependent on the training data, we decided to generate multiple gene predictions and build consensus gene models using EVidenceModeler (EVM) [START_REF] Haas | Automated eukaryotic gene structure annotation using EVidenceModeler and the program to assemble spliced alignments[END_REF] . The following gene prediction approaches were used: A) We trained Augustus with only the RNA-seq evidence and predicted gene models by taking into account RNA-seq, protein, Mikado and PASA assembled transcripts. B) We used the previously trained (in A) species model but with a modified weighting file (extrinsic.cfg) to give more weight to the protein evidence. C) We trained Augustus using both protein and RNA-seq evidence within Braker2 (EPT mode of Braker2). D) We used the RNA-seq evidence to automatically train genemark and obtain gene predictions. E) Finally, we run Augustus in the comparative mode with RNA-seq, transcript, and protein evidence and genome alignments inferred by mugsy [START_REF] Angiuoli | fast multiple alignment of closely related whole genomes[END_REF] . Generating this series of genome predictions was necessary as our preliminary analyses suggested that none of the predictions was superior but rather complementary. The proteomes used can be found in Supplementary Table 14.

Generating consensus gene models

We used EVM to derive consensus gene models best supported by the various evidence. We used all the previously generated gene predictions (gff files) and selected the best consensus gene models using protein (exonerate-mapped proteomes of species and the uniprot_sport plant dataset) and transcript evidence (Mikado and PASA assembled transcripts). We gave equal weights to each ab initio predictions, transcript evidence (weight 1), but increased the weight for Mikado loci (2) and PASA assembled transcripts (10). After deriving the consensus gene models, we used PASA and the PASA assembled transcripts to correct erroneous gene models, add UTRs, and predict alternative splice variants in two rounds.

Collinearity of the three hornwort genomes and collinearity across viridiplantae

We used the D-GENIES dot-plot tool [START_REF] Cabanettes | dot plot large genomes in an interactive, efficient and simple way[END_REF] with the default options to visually assess collinearity of the three genome assemblies. We also aligned the genomic sequences using the nucmer module of mummer [START_REF] Marçais | MUMmer4: A fast and versatile genome alignment system[END_REF] and assessed their differences using Assemblytics [START_REF] Nattestad | Assemblytics: a web analytics tool for the detection of variants from an assembly[END_REF] .

To study the collinearity across all plants, we first created orthogroups with 19 species' proteomes using Orthofinder2 [START_REF] Emms | OrthoFinder2: fast and accurate phylogenomic orthology analysis from gene sequences[END_REF] . The dataset included representatives from each major groups of land plants (Supplementary Table 14), and species experienced different numbers of large-scale duplication events 100 . Gff files and proteomes were retrieved from Phytozome v12 101 . We used I-ADHore3 102 to detect highly degenerate collinear blocks among bryophytes and vascular plants requiring a minimum of 3, 4 and 5 anchor points within each collinear region ( gap_size=30, cluster_gap=35, q_value=0.75, prob_cutoff=0.01, anchor_points=5, alignment_method=gg2, level_2_only=false) .

Identification of tandem repeats and centromeres

We run Tandem Repeats Finder 103 to identify tandem repeats with a minimum alignment score of 50 and a maximum period size of 2,000 bp. We then plotted repeat unit size against tandem array size to look for bimodal distribution. To localize centromeric regions in the A. agrestis Bonn genome, we generated dot-plots between a short-read-only assembly and the final chromosome-scale assembly. Because centromeric repeats are difficult to assemble using short-reads we expected that they will be missing from the Illumina assembly but will be present in the chromosomal-scale assembly. We also generated a self dot-plot of the A. agrestis Bonn genome to search for regions that are highly similar across scaffolds and are repetitive. Finally, we used the output of Tandem Repeats Finder 103 to search for tandem arrays with a period length of minimum 10bp and with a minimum tandem array length of 30 repeat units. We plotted the location of these tandem arrays along the chromosomes to visually assessed their distribution.

Screening for whole genome duplication

We used a combination of synonymous divergence (Ks) and synteny analyses to look for evidence of whole genome duplication in the Anthoceros genomes. For each genome, we used the DupPipe pipeline to construct gene families and estimate the age of gene duplications 104 . We translated DNA sequences and identified reading frames by comparing the Genewise 105 alignment to the best-hit protein from a collection of proteins from 25 plant genomes from Phytozome 101 . For each analysis, we used protein-guided DNA alignments to align our nucleic acid sequences while maintaining reading frame. We then used single-linkage clustering to constructed gene families and estimated Ks divergence using PAML 106 with the F3X4 model for each node in the gene family phylogenies. Because the Anthoceros genomes contain large numbers of pentatricopeptide repeat genes (PPR), we also repeated the analysis with all the PPR genes removed. PPR genes were identified based on the Orthofinder results (see below).

For synteny analysis, we used MCscan's "jcvi.compara.catalog ortholog" 107 function to search for and visualize intra-genomic syntenic regions, and used the "jcvi.compara.synteny depth" function for calculating syntenic depths. For comparison, we also carried out the same analysis for P. patens v3.3 and M. polymorpha v3.0 genomes; the former is known to have two rounds of whole genome duplications while the latter has none [START_REF] Bowman | Insights into land plant evolution garnered from the Marchantia polymorpha genome[END_REF][START_REF] Lang | The Physcomitrella patens chromosome-scale assembly reveals moss genome structure and evolution[END_REF] .

Transcription factor annotation

TAPs were annotated using TAPscan, according to Wilhelmsson et al [START_REF] Wilhelmsson | Comprehensive genome-wide classification reveals that many plant-specific transcription factors evolved in streptophyte algae[END_REF] and compared with selected other organisms using the major protein of each gene model (".1" splice variant). TF annotations were further manually checked and adjusted for annotation errors or missing annotations.

Gene family classification

We used Orthofinder2 99 to classify gene families of 25 plant and algal complete genomes, including the three hornworts reported here (Supplementary Table 15) into orthogroups. Orthofinder was ran using the default setting, except the "msa" option was used. A total of 31,001 orthogroups were circumscribed. The detailed gene count and classification results can be found in Supplementary Table 15.

Phylogenetic reconstruction of KNOX, LEAFY, WOX, and YABBY

For KNOX, AagrBONN.evm.model.Sc2ySwM.368.1986.6 was used as a query to BLASTp search at NCBI on 13th Sept 2019. The search database was nr limited to records that include: A. thaliana , Oryza sativa (japonica cultivar-group), Phalaenopsis equestris , Amborella trichopoda , Ceratopteris richardii , Selaginella moellendorffii , M. polymorpha , P. patens , K. nitens , Ostreococcus tauri , C. reinhardtii . The search parameters were otherwise as default. The hit sequences were downloaded and combined with the Anthoceros KNOX sequences, then aligned with FFT-NS-2 in MAFFT v7.427 108 . The alignment was manually inspected in Mesquite v3.6 109 and well conserved 149 sites of 51 sequences were included. Phylogenetic analysis based on Maximum Likelihood (ML) was conducted in MEGA X 110 . The best-fitting model was chosen as LG+G+I using the FindBestProteinModel function. A total of 100 bootstrap replicates were performed to evaluate branch support. "ML Heuristic Method" was set to "Subtree-Pruning-Regrafting -Extensive (SPR level 5)", and "No. of Discrete Gamma Categories" set to 5.

For LEAFY, AagrOXF evm.model.utg000049l.76.4 was used as a query to BLASTp search at NCBI on 30th Aug 2019. The search database was nr limited to records that include: A. thaliana , O. sativa (japonica cultivar-group), P. equestris , A. trichopoda , P. radiata , P. armandii , P. abies , C. richardii , S. moellendorffii , M. polymorpha , and P. patens . The search parameters were otherwise as default. The hit sequences were downloaded and combined with the Anthoceros LEAFY sequence and AHJ90704.1, AHJ90706.1, AHJ90707.1 from Sayou et al. 111 , then aligned with FFT-NS-2 in MAFFT v7.427 108 . The alignment was manually inspected and processed as described above to include 194 conserved sites of 20 sequences. Phylogenetic inference was done similarly as above but with LG selected as the best-fitting model.

For WOX, WOX genes in Anthoceros genomes were searched using the corresponding A. thaliana , P. patens and M. polymorpha proteins. Based on comparison among the three genomes, three gene models with excess intron predictions were manually revised, and one model was added. AagrOXF_evm.model.utg000018l.552.1 was used as a query to BLASTp search at NCBI on 9 October 2019. The search database was nr limited to records that include: A. thaliana , O. sativa (japonica cultivar-group), P. equestris , A. trichopoda , C. richardii , S. moellendorffii , M. polymorpha , P. patens , K. nitens , and C. braunii . Max target was set to 250 and the word size as 2. The search parameters were otherwise as default. The hit sequences were downloaded and combined with the Anthoceros WOX sequences, then aligned with einsi --maxiterate 1000 in MAFFT v7.429 108 . The alignment was manually inspected with Mesquite v3.6 109 and a matrix consisting of 58 included sites of 142 sequences were constructed. Sequences identical in the included region were treated as a single OTU during the phylogenetic analysis. The best-fitting model was chosen as JTT with ProteinModelSelection8.pl. The ML tree was inferred by RAxM 112 with -f a -\# 100 -m PROTGAMMAJTT and supplying -p and -x from random number generator. Bootstrap samples were generated with seqboot from PHYLIP package v3.697 113 and RAxML was run for each of them.

For YABBY , the 107 OTU dataset from Finet et al 114 was downloaded from treebase and combined with YABBY genes from Huperzia and Anthoceros . The sequences were aligned using einsi of MAFFT v7.450 108 . The aligned sequences were manually inspected with Mesquite and short sequences were removed and ambiguously aligned or gap containing sites were excluded. The best-fitting model was chosen as HIVB by ProteinModelSelection8.pl, and ML tree search followed what was described for WOX.

Phylogenetic reconstruction of stomata-related genes

An Anthoceros ICE/SCRM homolog sequence AagrBONN_evm.model.Sc2ySwM_368.1570.1 was used as a query to BLASTp search at NCBI on 7 Oct 2019. The search database was nr limited to records that include: A. thaliana , O. sativa (japonica cultivar-group), P. equestris , A. trichopoda , S. moellendorffii , P. patens , M. polymorpha , C. braunii , K. nitens . The word size was set to 2 and max target sequences as 250. The search parameters were otherwise set as the default. The hit sequences (100) were downloaded and combined with the Anthoceros ICE/SCRM sequences, then aligned with einsi --maxiterate 1000 in MAFFT v7.429 108 . The alignment was manually inspected with MacClade 4.08 and well-conserved 123 sites were included to result in alignment of 66 sequences. The sequence identical in the included region were treated as a single OTU during the phylogenetic analysis. The best-fitting model was chosen as JTTDCMUTF with ProteinModelSelection8.pl. The ML tree was inferred by RAxML with -f a -# 100 -m PROTGAMMAJTTDCMUTF and supplying -p and -x from random number generator. Bootstrap samples (1000 replicates) were generated with seqboot from PHYLIP package v3.697 113 and RAxML 112 was run for each of them. For ERECTA and TMM, the sequences of AagrOXF_evm.model.utg000083l.351.1 and AagrOXF_evm.model.utg000012l.100.1 were respectively used as the query and processed as in ICE/SCRM. Phylogenetic analyses were performed same as ICE/SCRM case, but with LG selected as the best-fitting model. For EPF and EPF-like gene family, we used the matrix compiled by Takata et al 115 and added the Anthoceros and M. polymorpha homologs. ML tree inference was done by IQ-TREE v1.6.1 with 1,000 replicates of UltraFast Bootstraps 116 .

Identification of orthologs to AMF symbiosis genes

Homologs to symbiotic genes were retrieved in 31 species covering the different plant lineages (Supplementary Table 11) using protein from the model plant Medicago truncatula and the tBLASTn v2.9.0+ 117 with a threshold e-value of 1e -10 . Sequences were aligned using MAFFT v7.407 108 with default parameters and alignments were cleaned using TrimAl v1.4 118 to remove positions with more than 20% of gaps. Resulting alignments were subjected to ML tree inference using IQ-TREE v1.6.1 119 . Prior to ML analysis, the best-fitting evolutionary model was tested using ModelFinder 120 and according to the Bayesian Information Criteria. Branch support was tested using 10,000 replicates of UltraFast Bootstraps 116 . Trees were visualized with the iTOL platform v4.4.2 121 .

Phylogenetic reconstruction of LCIB

The orthogroup OG0009668 was identified as the LCIB gene family containing C. reinhardtii LCIB -E genes. Additional hornwort LCIB homologs were retrieved from the 1KP transcriptome database 5 . To find other LCIB homologs, we ran BLASTp against the Phytozome database using both the Anthoceros and C. reinhardtii sequences as the query, and no hit could be obtained. Gene phylogeny was reconstructed based on the amino acid alignment done by MUSCLE 122 . IQ-TREE v1.6.1 119 was used to obtain the ML tree as outlined above.

F.-W.L., P.S., K.S., T.N. coordinated the project, M.S. carried out chromosome work, F.-W.L., P.S., T.N., D.H., S.C., G.K.-S.W. sequenced the genomes, F.-W.L., P.S. assembled the genomes, P.S., T.N. annotated the genomes, T.R., P.S. assembled and annotated organellar genomes, P.S., F.-W.L. performed synteny analyses, Z.L., M.S.B. performed Ks analyses, A.N., P.S. conducted RNA-seq experiment on developmental stages, F.-W.L., J.C.M. conducted RNA-seq experiment on cyanobacterial symbiosis, F.-W.L. conducted RNA-seq experiment on CO2 response, F.-W.L., M.W., A.K., I.D., P.S. analyzed RNA-seq data, N.P., S.R., M.W., K.S., P.S., E.F. characterized transcription factors, F.-W.L., P.S. performed gene family classification, J.K., P.-M.D. conducted analysis on AM symbiosis genes, I.M. conducted analysis on jasmonates, S.M., D.W. conducted analysis on auxin signaling, A.C. conducted analysis on ABA signaling, T.B. conducted analysis on strigolactone signaling, J.H.-G., M.A.B. conducted analysis on gibberellin signaling, S.dV. conducted analysis on salicylic acid signaling, J.dV., E.F. conducted analysis on genes associated with polyplastidy, J.H. conducted analysis on PIN proteins, S.W. conducted analysis on plastid targeted genes, E.F., T.N., F.-W.L. conducted analysis on stomatal development genes, F.-W.L., P.S., K.S., T.N., J.C.V., E.F., M.W. synthesized and wrote the manuscript.

Figure 1 .

 1 Figure 1. Genome organizations in the Anthoceros genomes . a , Circos plot of A. agrestis Bonn showing the densities of repeats, genes, and SNP with the A. agrestis Oxford genomes. No self-synteny can be found in the three Anthoceros genomes ( b-d ) nor in Marchantia polymorpha ( e ). Physcomitrella patens ( f ), on the other hand, shows a clear 1:1 and some 1:2 syntenic relationship, suggesting paleopolyploidy. In each panel the bar graph shows the proportion of the genome at different self-self syntenic levels, with the dotplot on the right.

Figure 2 .

 2 Figure 2. Transcription associated proteins (TAPs) and sporophyte development . a, The Anthoceros genomes have the smallest TAP repertoire among land plants. b, Gene expression profiles across different developmental stages in A. agrestis Bonn. Only selected genes known to regulate sporophyte and stomata development are shown here. c, Stomata of A. agrestis Bonn. d, Sporophytes (red arrowhead) and gametophytes (blue arrow) of A. agrestis Bonn.

Figure 3 .

 3 Figure 3. Evolution and genetics of symbiosis in hornworts . a , Orthologs of AMF symbiosis pathway genes can be found in hornworts, indicating their presence in the common ancestor of land plants. * Marchantia paleacea transcriptome was searched instead of M. polymorpha genome because the latter secondarily lost AMF. b , RAM1 is up-regulared during nitrogen starvation in both A. agrestis and A. punctatus . c , Reconstituted Anthoceros -cyanobacteria symbiosis. Arrowhead points to a cyanobacteria colony. d , Transcriptomic responses to nitrogen starvation and cyanobacterial symbiosis in A. agrestis . e , A suite of genes were highly up-regulated (>16 folds) under symbiosis in both A. agrestis and A. punctatus .

Figure 4 .

 4 Figure 4. Relationship between LCIB and pyrenoid-based carbon-concentrating mechanism . a , Phylogeny of LCIB . Numbers above branches are bootstrap support values (branches thickened when bootstrap > 70). b , Hornwort LCIBs have the conserved amino acid residues at the active site. Yellow and green arrowheads point to the zinc-binding and catalytic residues, respectively.

Table 1 .

 1 Assembly statistics of the three hornwort genomes.

		Estimated	Assembled	Contig/Scaffold	Contig/Scaffold	Assembly
		genome	genome	number	N50 length	approach
		size	size
	Anthoceros	124.5 Mb	116.9 Mb	1577/322	155.5Kb/17.3M	Illumina +
	agrestis				b/	Nanopore +
	Bonn					Hi-C
	Anthoceros	126.1 Mb	122.9 Mb	153/-	1.8Mb/-	Nanopore +
	agrestis					Illumina
	Oxford
	Anthoceros	129.4 Mb	132.8 Mb	202/-	1.7Mb/-	Nanopore +
	punctatus					Illumina

Acknowledgements

This project was supported by National Science Foundation DEB1831428 to F.-W.L., Swiss National Science Foundation (grants 160004 and 131726 ) to P.SZ., the Georges and Antoine Claraz Foundation (Switzerland) to P.SZ., The Forschungskredit and the University Research Priority Program "Evolution in Action" of the University of Zurich to M.W and P.SZ., IOS-1339156 and EF-1550838 to M.S.B., National Institute for Basic Biology (NIBB) Collaborative Research Program (13-710) to T.N., Japan Society for the Promotion of Science KAKENHI 26650143, 18K06367 to K.S., Short Term Postdoctoral Fellowship (PE14780) to E.F., Bill & Melinda Gates Foundation grant (OPP11772165) to P.-M.D., Spanish Ministry of Science, Innovation and Universities (BFU2016-80621-P) to J.H.-G. and M.A.B., European Research Council starting grant ("TerreStriAL") to J.dV., Foundation of German Business (sdw), Georges and Antoine Claraz Foundation, URPP Evolution in Action to A.N., German Science Foundation (WI4507/3-1) to S.W., Special Grant for Innovation in Research Program of the Technical University of Dresden (Germany) to D.Q. and S.W., Netherlands Organization for Scientific Research VICI grant (865.14.001) to D.W. and S.K.M. We thank Dr. Katsushi Yamaguchi and Prof. Shuji Shigenobu of Functional Genomics Facility at NIBB, Japan and Lucy Poveda, Catharine Aquino, Andrea Patrignani of the Functional Genomics Center Zurich (FGCZ) for sequencing support. Computational resources were partly provided by the Data Integration and Analysis Facility, NIBB and the NIG supercomputer at ROIS National Institute of Genetics.

Author Contributions


Single-molecule magnet behaviour in a Dy(III) pentagonal bipyramidal complex with a quasi linear Cl-Dy-Cl sequence. We report the synthesis and magnetic investigations of a dysprosium pentagonal bipyramidal complex [Dy(THF) 5 Cl 2 ][BPh 4 ] (1) exhibiting a linear Cl-Dy-Cl sequence suitable to provide a coordination environment allowing a zero-field slow relaxation of the magnetization. Besides, the complex shows also a dual luminescence originating from [BPh 4 ]  and Dy 3+ .

The ability for a coordination complex to present bistable magnetic states continues to catch attention of chemists and physicists since the discovery of this phenomenon nearly thirty years ago. [1][2] In this sense, the field of Single-Molecule Magnets (SMM) has taken an upturn with the use of lanthanide ions to design highly magnetically anisotropic complexes exhibiting exceptional performances, [3][START_REF] Tang | Lanthanide Single Molecule Magnets[END_REF][START_REF] Layfield | Lanthanides and Actinides in Molecular Magnetism[END_REF][START_REF] Ungur | [END_REF] making such molecular systems as potential candidates for applications in data storage, spintronics or quantum computing. 3, 7-9 In lanthanide SMM, the occurrence of a slow relaxation of the magnetization is directly correlated to the lanthanide crystal-field generated by the surrounding ligands that favours the creation of an anisotropic barrier, , opposing two antiparallel directions of the magnetic moment (±m J ). Intuitive synthetic strategies based on simple electrostatic considerations have been therefore applied to exploit the intrinsic magnetic anisotropy of 4f ions and engineer complexes displaying high energy barriers and blocking temperatures. Hence, one of the basic principles consists in stabilizing at most the electronic density of the 4f orbitals with a defined distribution of the negatively charged ligands. 10 This implies that oblate electronic density is efficiently stabilized by axial crystal-field, while prolate one demands an equatorial crystal-field. Such approach simultaneously increases the crystal-field splitting and minimizes the Quantum Tunnelling of the Magnetization (QTM), which provides direct relaxation within the ground state doublet. In addition to the relaxation of the magnetization through the first or higher excited states, systems with optimized molecular vibrations have to be engineered in order to reduce the spin-phonon coupling that induces alternative pathways (Raman or direct process) that clearly demand to be further comprehended.

3-6, 11

One of the most relevant approaches to design performing Dy 3+ -based SMM consists in stabilizing the Dy 3+ oblate electronic density by two ligands with high electronic density in axial positions and eventually other weak donor ligands in an equatorial plane. 6, 12 Indeed major breakthroughs in this field have recently been achieved by using coordination and organometallic chemistry [13][14][15][16][17][18] that have pushed back the limits to reach a magnetic hysteresis up to 60 K in a dysprosium metallocene complex 19-20 . On the other hand, it was shown that engineering pseudolinear complexes using simple coordination chemistry ligands is far to be straightforward owing to the difficulty to accurately control the crystal-field and symmetry around the lanthanide site. In this sense, all the aforementioned examples [13][14][15][16][17][18][19][20] used either negatively charged oxygen or organometallic basedligands. Surprisingly, chloride ion has been only rarely used to design complexes showing a linear Cl-Dy-Cl sequence 21-22 although it is known as the third most electronegative element after fluorine and oxygen. With this in mind, we report herein the synthesis, structure and magnetic studies of a dysprosium(III) dichloride pentagonal bipyramidal complex exhibiting a linear Cl-Dy-Cl sequence. The complex exhibits a zero-field slow relaxation of the magnetization associated with a genuine SMM behaviour and a dual luminescence.

Our design strategy relies on the formation of a cationic [DyCl 2 ]

+ moiety in association with a bulky anion to favour its crystallization. Thus, the reaction between DyCl The molecular structure obtained by X-Ray diffraction on single-crystal indicates that 1 crystallizes in the triclinic space group with two crystallographically independent complexes (noted Dy1 and Dy2) within the asymmetric unit. The coordination environment of both Dy 3+ sites is composed of two coordinated chlorides located in a trans fashion, as well as of five THF ligands in an equatorial plane leading to a seven coordination number (Fig. 1). The magnetic properties of 1 were investigated in both, static and dynamic modes by using a SQUID MPMS-XL. At 300 K, the T value of 14.20 cm 3 .K.mol 1 is in good accordance with the value of 14.17 cm 3 .K.mol 1 ( 6 H 15/2 ) expected for a unique Dy 3+ ion using the free-ion approximation. Lowering the temperature results in a monotonous decrease of T due to the thermal depopulation of the m J levels before dropping for temperature lower than 12 K (Fig. S2). The field dependence of the magnetization at 1.8 K shows a value of 6.08 Nβ under à 70 kOe field and does not show any saturation confirming the presence of a magnetic anisotropy (Fig. S2).

Besides, a sigmoidal shape is present at low fields, which has been ascribed in other systems to the occurrence of the high anisotropy or magnetic blocking.

14, 17

Fig. 2. Frequency dependence of the in-phase (') and out-of phase (") susceptibilities for 1 under a zero dc-field.

The presence of a SMM behaviour was investigated by alternate current measurements. The frequency dependence of the ac susceptibility under a zero dc-field reveals the presence of a single frequency dependent maximum in the out-of-phase susceptibility (") component indicating a slow relaxation of the magnetization (Fig. 2). Cole-Cole plots (Fig. S3) shows well-defined semi-circles which can be fitted with a generalized Debye model giving  parameter values ranging from 0.009 to 0.266 indicating a narrow distribution of relaxation times (Table S2). The increase of "observed at low temperature on the temperature dependence of " (Fig. S4) points out the presence of a QTM, which is further confirmed by monitoring the temperature dependence of the relaxation time, , which becomes temperature independent at low temperature (Fig. 3). Hence, the ln vs. T 1 plot shows a clear deviation from the linearity at low temperature confirming the presence of other relaxation processes. Consequently, the overall data range could be modelled using the following equation:

 1 =  0 1 exp(/kT) + CT m +  1
QTM (Eq. 1). 25 The first term accounts for a thermally activated process, while the second and third ones stand for two-phonon Raman and QTM, respectively. The m value was alternatively fixed to 5 (best fit) or 9, which are the possible values observed for Kramers ions. 26 The best fit parameters gives  = 58 ± 1 cm 1 ,  0 = (4.1 ± 0.9)  10 10 s, C = 0.038 ± 0.008 s 1 .K 5 and  QTM = 1.68 ± 0.03 ms, indicating that the magnetization relaxes through a combination of these three processes. To go further, the significant QTM could be shortcut by applying dc-field. As a consequence, the field dependence of the relaxation time at 6 K reveals the usual shape, which can be modelled with the equation  1 = DH 4 T + B 1 /(1+B 2 H²) + K (Eq. 2), for which the first term accounts for the direct process (for Kramers-ion), the second one for the QTM, while the K constant accounts for the field-independent Raman and thermally activated processes. The best fit parameters could be found in Table S3, while the optimum field (field at which the relaxation time is the highest) is in the range 1000-2000 Oe. The frequency dependence of the ac susceptibilities under a 1000 Oe field (Fig. S6) validates the reduction of the QTM, while the fitting of Cole-Cole plots (Fig. S7, Table S4) still shows low values of the  parameter (< 0.080) indicating a narrow distribution of the relaxation times. Fitting of the temperature dependence of the relaxation time could be performed with the following model:

 1 =  0 1 exp(/kT) + CT m + AT (Eq. 3)
(Fig. 3) in which the third term accounts for the direct process leading to the following fit parameters:  = 65.2 ± 0.2 cm 1 ,  0 = (1.145 ± 0.002)  10 13 s, C = 0.004 ± 0.002 s  1.K 5 (the A parameter is found negligible). All these results point out a modest value of the anisotropic barrier, which could be rather surprising giving the quasilinearity of the Cl-Dy-Cl sequence. Evaluation of the orientation of the anisotropic axes from the MAGELLAN 27 software based on electrostatic considerations and for both independent complexes confirms that this latter is directly oriented along the Cl-Dy-Cl string with anisotropic axis-Dy-Cl angle of around 2.0 and 0.4° for Dy1 and Dy2 sites, respectively (Fig. S8). Although the negatively charged chloride anions in the trans position impose the orientation of the anisotropic axis from a pure electrostatic point of view, the presence of strongly interacting oxygen in the pentagonal plane provides a transverse magnetic anisotropy component and therefore a decrease of the overall anisotropy, which explains the modest value of the anisotropic barrier. Additionally, the relatively weak crystal-field generated by the chloride ions may not be able to afford a large m J splitting. However and by comparison, a benzoquinone(BQ)-bridged dinuclear dysprosium complex [BQDy 2 Cl 2 (THF) 6 ] showing a close environment with respect to 1 does not exhibit a slow relaxation of its magnetization due to the strong equatorial field imposed by the BQ ligands that induces a rotation of the magnetic axis in the basal plane. 22 On the other hand, a dysprosium complex based on a tridentate NCN pincer ligand [Dy(NCN)Cl 2 (THF) 2 ] with two chlorides in trans configuration shows an SMM behavior with a relaxation through the second excited state. In this example the magnetic axis is imposed by a strong Dy-C interaction within the equatorial plane.

21

The optical properties of 1 were investigated. The absorption spectrum of 1 in THF shows a broad band at 270 nm and ascribed to the *   transitions of the [B(Ph) 4 ]  moiety (Fig. S9). 28 Then, the photoluminescence properties of 1 were investigated in solid state. The room temperature emission spectra acquired at different excitation wavelengths show a dual luminescence constituted by an excitation dependent broad band in the 410-435 nm region and the typical Dy 3+ -based luminescence at 485, 577 and 668 nm. The former broad band most likely originates from the tetraphenylborate moiety. Noticeably, the luminescence of tetraphenylborate salts is known to be relatively complex since this arises from the formation of excitons. [29][30] The excitation spectrum (Fig. S10) monitored in the intra 4f In conclusion, a simple dysprosium pentagonal bipyramidal complex exhibiting a linear Cl-Dy-Cl arrangement reveals a slow relaxation of the magnetization associated with a genuine SMM behaviour and a dual luminescence making the compound a bifunctional system. Despite the linear arrangement, a relatively modest anisotropic barrier was observed, which can be directly related to the presence of strongly coordinated THF solvates in the pentagonal plane that decrease the anisotropy. Yet, substitution of such solvates with softer Lewis bases may greatly enhance the magnetic properties.

  3 and NaBPh 4 was carried out in THF at room temperature in anhydrous conditions and allowed to obtain single crystals of [DyCl 2 (THF) 5 ][B(C 6 H 5 ) 4 ]•0.5(THF) (1) (Scheme S1) after slow diffusion of hexane.

Fig. 1 . 1 .

 11 Fig. 1. Top: Molecular structure of one of the independent complex [DyCl2(THF)5] + in 1. Colour code: orange, Dy; red, O; grey, C; green, Cl. Hydrogen atoms have been omitted for clarity. Bottom: View showing the pentagonal bipyramidal arrangement.

Fig. 3 .

 3 Fig. 3. Temperature dependence of the relaxation time for 1 using the ac susceptibility data at 0 Oe and 1000 Oe. The solid line represents the fit with Eq. 1 (0 Oe) and Eq. 3 (1000 Oe).

Fig. 4 .

 4 Fig. 4. Room temperature solid-state emission spectra of 1 acquired at different excitation wavelengths. The intra 4f 9 Dy 3+ transitions are reported.

  In addition, a broad excitation band ascribed to the [B(Ph) 4 ]  entity could be also observed at 300 nm. Lowering the temperature (77 K) results in the increase of the relative of the Dy 3+ -based luminescence with respect to the broad emission band. These results indicate therefore than 1 exhibit a characteristic dual luminescence.

	-8 -7 -6 -5 -4 -3 Lns	0 Oe 1000 Oe	Intensity / arb. units	6 H 15/2	4 F	9/2	6 H	13/2	   	exc = 300 nm exc = 325 nm exc = 365 nm exc = 350 nm
	-9								4 F 9/2	6 H 11/2
		9 4 F 9/2  3+ transitions reveals the presence of sharp lines at 325, 350 6 H 13/2 Dy and 365 nm pointing a direct sensitization within the 4f transitions. 0.2 0.3 0.4 0.5 T -1 400 450 500 550 600 650 700 750 -1 / K Wavelength / nm

intensities

Acknowledgements

The financial support of the Russian Science Foundation is highly acknowledged (Project № 17-73-30036). The French authors thank the University of Montpellier, CNRS and PAC of ICGM.

Notes and references

Conflicts of interest

There are no conflicts to declare.