Sex Determination by Two Y-linked Genes in Garden Asparagus

The origin and early evolution of sex chromosomes has been hypothesized to involve the linkage of factors with antagonistic effects on male and female function. Garden asparagus ( Asparagus officinalis L.) is an ideal species to investigate this hypothesis, as the X and Y chromosomes are cytologically homomorphic and evolved from an ancestral autosome pair in association with a shift from hermaphroditism to dioecy. Mutagenesis screens paired with single-molecule fluorescence in situ hybridization (smFISH) directly implicate Y-specific genes that respectively suppress female (pistil) development and are necessary for male (anther) development. Comparison of contiguous X and Y chromosome assemblies shows that hemizygosity underlies the loss of recombination between the genes suppressing female organogenesis ( SUPPRESSOR OF FEMALE FUNCTION, SOFF ) and promoting male function ( TAPETAL DEVELOPMENT AND FUNCTION1, aspTDF1 ). We also experimentally demonstrate the function of aspTDF1. These findings provide direct evidence that sex chromosomes can function through linkage of two sex determination genes.


INTRODUCTION
Beginning with Nettie Stevens' (Stevens, 1905) elegant explanation of male and female mealworm gamete differences involving "sex chromosomes", sex chromosomes have been characterized in dioecious species (with separate male and female sexes) across all major eukaryotic lineages. Flowering plants offer a unique perspective to the origin and early evolution of sex chromosomes across the tree of life, given that dioecy has evolved hundreds of independent times from hermaphroditic ancestors across the angiosperm phylogeny (Renner, 2014).
Several models have been proposed to explain the transition from a hermaphroditic species with autosomal chromosomes to a dioecious plant species in association with the evolution of an X/Y or Z/W sex chromosome pair. One model developed over time by both Westergaard (Westergaard, 1958;Charlesworth, 2018) and Charlesworth and Charlesworth (Charlesworth and Charlesworth, 1978) hypothesizes that a sex chromosome could evolve from an autosomal pair via mutations in two perfectly linked sexually antagonistic genes in a non-recombining region of a proto-Y chromosome: one gene that promotes anther development, and another that suppresses female organ development. Westergaard supported his two-gene model with genetic data from papaya (Carica papaya) (Storey, 1953), and his own work in Silene (formerly Melandrium) (Westergaard, 1953(Westergaard, , 1958. However, in the genomics age it has been notoriously difficult to identify causal sex determination genes that fit this model in plants and animals, given the difficulty of assembling sex chromosome sequences, particularly in older sex chromosomes that have degenerated, accumulated repetitive DNA and heterochromatin, undergone structural variations, and possibly experienced gene turnovers (Charlesworth, 2019). Despite these complications, several sex determination genes have recently been identified in plants, including a non-coding RNA in persimmon (genus Diospyros) (Akagi et al., 2014(Akagi et al., , 2016(Akagi et al., , 2018, and two genes on the kiwifruit (genus Actinidia) Y chromosome (Akagi et al., 2019). Garden asparagus (Asparagus officinalis L.) has long been identified as a potential model system for investigating the genetic basis of sex determination in dioecious plant species (e.g. Irish and Nelson, 1989;Dellaporta and Calderon-Urrea, 1993;Westergaard, 1958). Here we extend our earlier genomic analyses of the Y-chromosome in garden asparagus (Harkess et al. 2017) to test the hypothesis that sex determination in garden asparagus is based on a two-gene, male-promoting and female-suppressing genetic M-locus (or Sex Determination Region) on the Y chromosome (Marks, 1973) and elucidate the mechanism for loss of recombination across the multi-gene Mlocus. New functional analyses, comparisons of the X and Y chromosomes, and single-molecule in situ data support two Y-specific genes as being sufficient for sex determination in garden asparagus.

RESULTS AND DISCUSSION
Two Y-linked genes determine sex in garden asparagus X and Y sex chromosomes and dioecy appear to have evolved from hermaphroditic ancestors once (Kubota et al., 2012) or potentially twice in the Asparagus genus (Norup et al., 2015). Sequencing of a YY "supermale" garden asparagus genome revealed the existence of a nearly one megabase non-recombining, male specific sex determination region (SDR) on the Y chromosome with only thirteen gene models (Harkess et al., 2017), two of which we identify here as the master sex determination genes ( Figure 1A). Previous work to identify the sex determination genes in the non-recombining regions of the garden asparagus Y chromosome utilized gamma irradiation and spontaneous mutants to identify a Y-specific, dominant female suppression gene, SUPPRESSOR OF FEMALE FUNCTION (SOFF) (Harkess et al., 2017). Gamma irradiation knockouts of the entire ~1-megabase (Mb) Y-linked sex determination region (Y-SDR) resulted in male-to-female conversion, whereas a single gene knockout of the Y-specific SOFF gene converted XY males to hermaphrodites with functioning styles and receptive stigmas ( Figure 1C). An additional spontaneous frameshift mutant verified SOFF as the female-suppressor. As predicted by the two-gene models for a dominant M-locus on the Y-chromosome (Westergaard, 1958;Charlesworth, 2018;Marks, 1973;Charlesworth and Charlesworth, 1978), a homolog of DEFECTIVE IN TAPETAL DEVELOPMENT AND FUNCTION1 (TDF1), encoding an R2R3 MYB transcription factor for which knockouts in Arabidopsis thaliana are male-sterile (Zhu et al., 2008), was identified in the non-recombining Y-specific region of the asparagus Y-chromosome (Harkess et al., 2017) and found to be male-specific across several, but not all, dioecious Asparagus species (Murase et al., 2017;Tsugama et al., 2017). The garden asparagus aspTDF1 was able to rescue the tdf1 mutation in Arabidopsis (Murase et al., 2017;Tsugama et al., 2017), but the function of aspTDF1 has not been validated in

Asparagus.
Here we exposed an all-XY male garden asparagus seed population to ethyl methanesulfonate (EMS) mutagenesis and recovered a male-to-neuter individual that lacked fully developed anthers ( Figure 1B).
Resequencing of the Y-specific aspTDF1 gene in this mutagenized individual revealed a single nucleotide mutation that induced a premature stop codon in the predicted aspTDF1 protein-coding sequence (Supplemental Figure 1). Thus we have now demonstrated that functional modification of one or both of the sex determination aspTDF1 and SOFF genes can convert XY males into one of three different sexual forms: knockout of the SOFF gene converts males to hermaphrodites, knockout of the aspTDF1 gene converts males to neuters, and knockout of both TDF1 and SOFF converts males to females ( Figure 1C). Knockouts for these two sexually antagonistic Y-linked genes, one that promotes anther formation (TDF1) and another that suppresses pistil development (SOFF), support the long-standing "two-gene" model of sex determination in dioecious species as posited by Westergaard (Westergaard, 1958) and Charlesworth and Charlesworth (Charlesworth and Charlesworth, 1978).

The structure of the X chromosome
One necessary feature of plant sex chromosomes is a mechanism to suppress recombination between the heterologous pair (e.g., between X and Y, or Z and W), maintaining the sex determination genes in linkage disequilibrium. Given the multiple independent origins of dioecy and sex chromosomes across flowering plants (Renner, 2014), it is unsurprising that suppressed recombination between X/Y or Z/W can initiate and persist via multiple mechanisms. For example, in papaya, the non-recombining region of the Y chromosome has repeatedly inverted and expanded across the already non-recombining centromere (Wang et al., 2012). In poplar trees (Populus sp.), a translocation event of a sex-specific region to the distal tip of chromosome 19 inhibits recombination between the sex chromosome pair (Yin et al., 2008;Tuskan et al., 2012). These initial events may be followed by extensive Y-linked repeat proliferation, exemplified by the large and expanded heterochromatic Y chromosomes in Silene latifolia and Coccinia grandifolia (Sousa et al., 2017).
Previous inferences about the structure and composition of the garden asparagus X chromosome have relied on resequencing data aligned to the YY reference genome (Harkess et al. 2017). We generated a de novo whole genome assembly of a doubled-haploid XX garden asparagus individual that is a sibling of the YY genotype for which we have a published chromosomal genome assembly (Harkess et al. 2017;Supplemental Figure 2). Assembly of roughly 32X coverage (estimated genome size 1.3 Gb) PacBio RSII reads yielded a 1.16-gigabase (Gb) genome assembly with 5,966 contigs (contig N50 = 385 kb), and further scaffolding with a BioNano Genomics Irys BspQI optical map reduced the scaffold count to 4,502 (scaffold N50 =1.67 Mb). The assembly was anchored against the YY genome pseudomolecules ( Figure 1D). BUSCO metrics (Simão et al., 2015) show that 85.3% of the full-length embryophyte database orthologs are present in the XX genome assembly, similar to the previous YY assembly (Supplemental Figure 3), suggesting that both assemblies share similar degrees of completeness.
The X chromosome assembly contained a contiguously assembled 163-kb region of X-specific sequence flanked by the same genes surrounding the nearly 1-Mb Y-specific sex determination region on the Y chromosome ( Figure 1E). The co-linearity of the flanking regions and absence of X/Y gametologs suggest that hemizygosity rather than an inversion or translocation is responsible for the lack of recombination between Xand Y-specific SDRs in garden asparagus. The SOFF and aspTDF1 genes are not present on the homologous X region or elsewhere in the XX female genome. There is no evidence of repeated inversions leading to the formation of strata on the asparagus Y chromosome ( Figure 1E), as has been inferred for older and more degenerate plant sex chromosomes like Carica papaya (Wang et al., 2012) and Silene latifolia (Bergero et al., 2008;Kazama et al., 2016). However, the finding that aspTDF1 is not sex-linked in all dioecious Asparagus species (Harkess et al. 2017, Murase et al., 2017Tsugama et al., 2017) indicates that gene content in the sexdetermination region, and indeed the molecular basis of sex determination, are evolving in the genus.
One possible consequence of sex chromosome evolution is that genes with male function will accumulate on the Y chromosome over time, and likewise genes with female-biased function will accumulate on the X chromosome. Translocations of entire sex-determination regions are even possible (Tennessen et al., 2018). The 163-kb X-linked region contains five expressed transcripts: three are annotated as uncharacterized proteins, one is a chloroplast-encoded gene shared with the Y chromosome, and one is X-linked with no Y homolog (Supplemental Table 1). This annotated X-linked gene is a homolog of Arabidopsis NO TRANSMITTING TRACT (WIP2/NTT), which encodes a CH2H/C2HC zinc finger transcription factor that is specifically expressed in the transmitting tract and funiculus of ovules, playing a key role in the development of fruit (Marsch-Martínez et al., 2014;Crawford et al., 2007).
Knockouts of WIP2/NTT in A. thaliana result in plants that inhibit the movement of the pollen tubes through the carpel and into ovules, severely reducing or eliminating male fertility. The development of functional ovules and seeds in YY males with mutated or suppressed SOFF genes suggests that the X-linked aspWIP2/NTT is not necessary for female function in greenhouse-grown plants, but controlled experiments have not been done to test the relative seedset of asparagus females with and without intact X-linked aspWIP2/NTT genes.

Single molecule FISH expression quantification
Little is known about the mode of action of the SOFF protein, other than it contains a PFAM annotation "Domain of Unknown Function 247 (DUF247)". Previous attempts at RNA-seq of whole spear tips (comprised of a mixture of both vegetative tissue and floral buds) were able to clearly detect the expression of aspTDF1, but not SOFF (Harkess et al., 2015). Intriguingly, a DUF247 domain-containing protein co-segregates with the S Self-Incompatibility locus in perennial ryegrass (Lolium perenne) (Manzanares et al., 2016), hinting at the possibility of a conserved potential suppressive function.
To clarify the spatial expression patterns of both genes, and to better elucidate their respective functions, single-molecule fluorescence in situ hybridization (sm-FISH) was used to track the expression and localization of the sex-determination genes on XX and YY spear tips. Since neither of the two sex determination genes is present in females, the XX spear tips serve as a negative control. EF1-alpha expression is used as a positive control in both XX and YY, and indeed is expressed ubiquitously throughout all sections. As expected from Arabidopsis (Zhu et al., 2008), in situs show that aspTDF1 expression is limited to the anther tapetal layer in developing YY supermale buds, and that as expected there is no detectable expression in XX female buds ( Figure 2A). SOFF is weakly expressed in YY supermale buds, and like aspTDF1, is not present in XX females ( Figure 2A). To resolve the expression levels of SOFF and aspTDF1 in flower buds, the smFISH signals for both aspTDF1 and SOFF were quantified through a Z-stack of entire developing flower buds. The mean expression of aspTDF1 (288 copies) and SOFF (99 copies) through a Z-stack of single flower buds highlights the efficiency of SOFF in inhibiting female organogenesis even in low mRNA copy number relative to aspTDF1 ( Figure 2B; Supplemental Table 2). Overall, while there is clear tapetum-specific expression of aspTDF1, SOFF shows less clear spatial gene expression patterns, and requires a finer staging of developing gametophyte tissue to more deeply explore its function.

Conclusions
Dioecious garden asparagus has an X/Y sex chromosome pair that determines sex via the action of two sexually antagonistic genes in linkage disequilibrium on the Y chromosome: a suppressor of female organogenesis (SOFF) and a promoter of anther development (aspTDF1). Though long hypothesized in several species, we show functional evidence that the evolution of two sex determination genes on a Y chromosome can determine plant sex, supporting the long-standing hypothesis posited by Mogens Westergaard (Westergaard, 1958) and Brian and Deborah Charlesworth (Charlesworth and Charlesworth, 1978). Earlier work (Harkess et al. 2017;Murase et al., 2017;Tsugama et al., 2017) indicated that aspTDF1 was not Y-specific in the last common ancestor of all dioecious Asparagus species. Together these findings suggest that whereas aspTDF1 in necessary for male function, either there was another ancestral male promoter gene responsible for the origin of dioecy in the genus, or there were multiple origins of dioecy through linkage of alternative male promoting genes to SOFF in dioecious lineages with autosomal aspTDF1 genes. Assessment of these hypotheses will require comparative analyses of Y chromosomes in multiple dioecious Asparagus species and orthologous autosomal chromosomes in closely related hermaphrodite species.

EMS mutagenesis
Twelve hundred seeds of the all-male hybrid cultivar Marte (courtesy of Dr. Agostino Falavigna, Agricultural Research Council, Lodi, Italy) were treated using 0.4 % (v/v) ethyl methanesulfonate (EMS). Seeds were imbibed for 6 h in rolling flasks filled with 100 mM Na2HPO4. This buffer was replaced by Na2HPO4 buffer plus 0.4% EMS, and seeds were further incubated in rolling flasks for 15 h at room temperature. EMS treatment was stopped by washing the seeds two times for 20 min each in 100 mM Na2S2O3 neutralizing buffer.
After neutralization, seeds were rinsed three times in demineralized water and incubated for another 10 min in rolling flasks with demineralized water. Subsequently, the seeds were dried overnight on filter paper in a laminar airflow cabinet at half air speed.
Plants were grown in the field at Limgroup in early Spring, where a neuter shoot was discovered on a M1 plant by visual inspection. Meristems of both the WT and mutant shoots taken from the same plant were grown in vitro and rooted as previously described (Qiao and Falavigna, 1989) to obtain ten cloned plants of each, all of which retained the original phenotype. Sanger sequencing of the TDF1 gene was performed with primer pairs (Forward: AGCCCTTAAGGTTAAATGTCG; Reverse: CATCATGATATAAAATCCTCAATCAAA).

Asparagus officinalis XX female genome assembly and annotation
High molecular weight DNA was isolated from a doubled haploid XX female sibling to the sequenced YY asparagus genome (Harkess et al., 2017) using the Bionano Genomics IrysPrep Plant kit. A 40-kb PacBio library was constructed and size-selected on a SageScience BluePippin, and then sequenced with P6C4 chemistry and six-hour movie lengths on a Pacific Biosystems RSII. A total of 45 SMRT cells were sequenced.
Raw data are deposited in NCBI BioProject PRJNA603233. Subreads greater than 8 kb in length were assembled with FALCON (v0.3.0). The same high molecular weight DNA was used to generate a BioNano Genomics optical map labeled with Nt.BspQI and imaged to 70X coverage. The assembly was scaffolded against the optical map using the HybridScaffold program in Irysview with default options. The scaffolded assembly was polished using nearly 40X coverage Illumina pair-end 100 nt reads from the same individual with Pilon v1.22 (Walker et al., 2014), and then assembled into pseudomolecules by comparison to the YY reference genome using CoGE Syntenic Path Assembly (Lyons, 2008). The assembly quality was estimated using BUSCO v3.1.0 (Simão et al., 2015) with the embryophyte odb9 database, trained against the maize annotation set.
The assembly was annotated for expressed transcripts using RNA-seq reads from spear tip tissue (which includes developing floral buds) of four XX females (NCBI BioProject 259909). Reads were aligned to the genome with STAR v2.5.3a (Dobin et al., 2013) using default options and "--outSAMstrandField intronMotif", and transcripts were assembled using Stringtie v1.3.3b (Pertea et al., 2015) with default options. A single Xlinked Pacbio contig spanned the two pseudoautosomal (PAR) boundaries of the Y. Sequence comparison of the X and Y chromosomes was performed using MUMMER v3 Nucmer (Delcher et al., 2003).

Sample preparation for in situ
Asparagus flower buds of various sizes were cut from the spear and fixed in a 20 ml glass vial using 4% paraformaldehyde in 1x PHEM buffer (5 mM HEPES, 60 mM PIPES, 10 mM EGTA, 2 mM MgCl2 pH 7).
Samples were fixed three times (15 min each) in a vacuum chamber at 0.08 MPa. After fixation, samples were sent for paraffin embedding at the histology lab in the Nemours/Alfred I. duPont Hospital for Children (Wilmington, DE).

Fluorescence in situ hybridization
Garden asparagus bud samples were sectioned using a paraffin microtome and dried on poly-l-lysine coated 22 x 22 mm #1.5 coverslips (Carl Zeiss Microscopy, LLC, Cat# 474030-9020-000). Samples were then ) and re-hydrated by going through an ethanol series of 95, 80, 70, 50, 30, 10% (vol/vol) (30 s each) and water (1 min) at room temperature. After protease (Sigma, P5147) digestion (20 min, 37°C), samples were treated with 0.2% glycine (Sigma-Aldrich, G8898) for 2 min, followed by a TEA treatment with 1.3% (v/v) triethanolamine (Sigma-Aldrich, 90279), 0.4% (v/v) HCl, and 0.5% (v/v) acetic anhydride (Sigma-Aldrich, A6404). After two washes in 1x PBS buffer (phosphate-buffered saline), samples were dehydrated and then hybridized with smFISH probes. The smFISH probes were designed to bind specifically across the length of the target RNAs (Markey et al., 2014). Briefly, 40 probes (each 20 nt long) were designed for each target, and each probe was synthesized with 3' amino modification (LGC Biosearch Technologies, CA). All the probes of one set were pooled and en mass coupled with TMR (tetramethylrhodamine) or Texas Red and the labeled probe fraction was purified using HPLC (Batish et al., 2011). The probes were diluted in hybridization buffer containing 10% formamide, yeast tRNA, dextran sulfate and RNAse inhibitor to make hybridization mix. The samples were hybridized with hybridization mix overnight in a humid chamber at 37°C. The samples were washed two times with 1X SSC buffer (sodium saline citrate) containing 10% formamide and a final wash was done using 1x TBS buffer (Trisbuffered saline). Samples were mounted using SlowFade™ Diamond Antifade Mountant (ThermoFisher Scientific, S36967).
Fluorescence was detected by spectral unmixing of autofluorescence spectra using laser scanning confocal microscopy on a Zeiss LSM 880 multiphoton confocal microscope. A 15% 561 nm laser and an alpha Plan-Apochromat 100x/1.46 oil lens were used to acquire TDF and SOFF mRNA (with probes labeled with Texas Red) in situ images using the online figure print mode. Pure TxRed was used as positive control, and autofluorescence from asparagus tissue was used as negative control for spectral bleed-through. The 555 nm laser was used to detect positive control EF1a mRNA labeled with TMR-linked probes, with the alpha Plan-Apochromat 100x/1.46 oil lens, in spectra mode. After images were taken, each image was spectra unmixed.
Spectral data for the pure TMR fluorophore were used as positive control, and non-labeled samples were used as negative control spectra for spectral bleed-through. The brightness and contrast of images in the same figure panel were adjusted equally and linearly in Zen 2010 (Carl Zeiss).

Image quantification
Localization events were counted using Volocity (Perkin Elmer). For each image, noise was first removed using the default value from Volocity, and then the localization events were quantified using a pipeline consist of Find Object (in the channel with in situ signal), Separate Touching Objects (0.05 µm 3 ), Exclude object (<0.05 µm 3 ), and Exclude Object (>0.8 µm 3 ). Five replicates for each Z-stack of images were used for calculating the copy number for each gene. Z-stacks represent stacks of images through a focal plane. A p-value was calculated using a t-test assuming equal variance.

Accession Numbers
The Asparagus officinalis XX female genome and browser is available at CoGE (id35894; https://genomevolution.org/coge/GenomeView.pl?embed=&gid=35894). Raw PacBio RSII data is available at BioProject PRJNA603233.  (Harkess et al., 2017) and this manuscript's XX genome assemblies. Table 1: X-linked BLASTX (nr) gene annotations. Table 2: T-test statistics. Gamma irradiations that delete the entire Y-SDR region convert XY males into females. Single-gene gamma irradiation mutants of SOFF convert XY males into hermaphrodites. A premature stop codon in aspTDF1 converts XY males into neuters. (D) Synteny between the XX PacBio plus BioNano genome assembly scaffolds and the 10 assembled chromosomes for the published YY genome (Harkess et al., 2017). (E) Microsynteny of the X chromosome against the Y chromosome (blue lines) with the non-recombining regions of both chromosomes (red blocks). Blue lines indicate MUMMER alignment matches greater than 1.5 kb with 90% minimum identities.