Author response: A unicellular relative of animals generates a layer of polarized cells by actomyosin-dependent cellularization

Abstract

Article Figures and data Abstract Introduction Results Discussion Materials and methods Data availability References Decision letter Author response Article and author information Metrics Abstract In animals, cellularization of a coenocyte is a specialized form of cytokinesis that results in the formation of a polarized epithelium during early embryonic development. It is characterized by coordinated assembly of an actomyosin network, which drives inward membrane invaginations. However, whether coordinated cellularization driven by membrane invagination exists outside animals is not known. To that end, we investigate cellularization in the ichthyosporean Sphaeroforma arctica, a close unicellular relative of animals. We show that the process of cellularization involves coordinated inward plasma membrane invaginations dependent on an actomyosin network and reveal the temporal order of its assembly. This leads to the formation of a polarized layer of cells resembling an epithelium. We show that this stage is associated with tightly regulated transcriptional activation of genes involved in cell adhesion. Hereby we demonstrate the presence of a self-organized, clonally-generated, polarized layer of cells in a unicellular relative of animals. Introduction Cellularization of a coenocyte — a multinucleate cell that forms through sequential nuclear divisions without accompanying cytokinesis — is a specialized form of coordinated cytokinesis that results in cleavage into individual cells. Cellularization commonly occurs during development of animals, plants and unicellular eukaryotes (Hehenberger et al., 2012; Mazumdar and Mazumdar, 2002). Despite the similarities, distinct mechanisms are involved in the cellularization of these coenocytes. During endosperm development in most flowering plants, coenocytes cellularize through cell wall formation around individual nuclei, forming a peripheral layer of cells surrounding a central vacuole (Hehenberger et al., 2012). This is coordinated by the radial microtubule system and is dependent on several microtubule-associated proteins (Pignocchi et al., 2009; Sørensen et al., 2002). In apicomplexan parasites, cytokinesis of multinucleate schizonts occurs by budding, which is driven by a polarized microtubule scaffold (Francia and Striepen, 2014). In contrast, in a model animal coenocyte, the syncytial blastoderm of the fruit fly Drosophila melanogaster, cellularization is accomplished through plasma membrane invaginations around equally spaced, cortically positioned nuclei (Farrell and O'Farrell, 2014; Mazumdar and Mazumdar, 2002). This process relies on extensive membrane remodeling (Lecuit, 2004; Sokac and Wieschaus, 2008; Figard et al., 2016) regulated by zygotically transcribed genes (Schweisguth et al., 1990; Schejter and Wieschaus, 1993; Postner and Wieschaus, 1994; Hunter and Wieschaus, 2000; Lecuit et al., 2002) and driven by microtubules and a contractile actomyosin network (Mazumdar and Mazumdar, 2002). It depends on several actin nucleators, such as the Arp2/3 complex (Stevenson et al., 2002) and formins (Afshar et al., 2000). It also requires multiple actin-binding proteins, including myosin II (Royou et al., 2002), which mediates actin cross-linking and contractility, as well as septins (Adam et al., 2000; Cooper and Kiehart, 1996), cofilin (Gunsalus et al., 1995) and profilin (Giansanti et al., 1998). In addition, it depends on cell-cell adhesion proteins including cadherin, and alpha- and beta-catenin (Hunter and Wieschaus, 2000; Wang et al., 2004). This coordinated cellularization results in the formation of a single layer of polarized epithelial tissue, also known as cellular blastoderm (Mazumdar and Mazumdar, 2002). This actomyosin-dependent cellularization is common in early insect embryos and also commonly observed in germline development of many animals (Haglund et al., 2011), such as the nematode C. elegans (Priti et al., 2018), however, whether this mechanism of cellularization is found outside animals, remains unknown. Among holozoans — a clade that includes animals and their closest unicellular relatives (Figure 1A) — ichthyosporeans are the only known lineage besides animals that forms coenocytes during their life cycles (Mendoza et al., 2002; de Mendoza et al., 2015). All characterized ichthyosporeans proliferate through rounds of nuclear divisions within a cell-walled coenocyte, followed by release of newborn cells (Marshall and Berbee, 2011; Marshall and Berbee, 2013; Suga and Ruiz-Trillo, 2013; Whisler, 1968). We have previously suggested that they undergo cellularization (Ondracka et al., 2018; Suga and Ruiz-Trillo, 2013). However, at present nothing is known about the ichthyosporean cellularization, and whether it involves animal-like mechanisms. Figure 1 with 1 supplement see all Download asset Open asset Cellularization dynamics in Sphaeroforma arctica. (A) Phylogenetic position of the ichthyosporean Sphaeroforma arctica in the tree of life. (B) Time-lapse images of the life-cycle of S. arctica show cell-size increase prior to cellularization and release of new-born cells. Associated with Video 2. Bar, 50 µm. (C) Kymographs of 4 distinct cells undergoing cellularization with the time of flip (yellow) and cell release (orange) indicated for each. An example of the advancing dark front is shown for cell 4. Bar, 50 µm. (D) Mean cell diameter over time of 65 cell traces aligned to Flip reveals distinct cell stages: Growth, Pre and post-flip and cell release. (E) Duration of growth, cellularization and post-flip represented as box-plots (N°cells > 100). MAD (Median absolute deviation) over median is used as a measure of variability. Here, we comprehensively characterized cellularization in the ichthyosporean Sphaeroforma arctica, in which we have previously shown that coenocytic cycles can be synchronized (Ondracka et al., 2018). We used imaging, transcriptional profiling, and pharmacological inhibition, to study the gene expression dynamics, morphological rearrangements, and mechanisms of cellularization. We found that cellularization is accomplished by inward plasma membrane invaginations driven by an actomyosin network, forming a polarized layer of cells. Time-resolved transcriptomics revealed sharply regulated expression of cell polarity and cell adhesion genes during this stage. Finally, we show that this process depends on actin nucleators and Myosin II, and we reveal the temporal order of the actomyosin network assembly. Together, these findings establish that cellularization of ichthyosporeans proceeds by mechanism conserved between animals and ichthyosporeans. Results Growth and cellularization in S. arctica are temporally separated stages of the coenocytic cycle To determine the timing of cellularization in synchronized cultures, we established long-term live imaging conditions. Individual coenocytes exhibit growth in cell size until approximately 60 hr, after which they undergo release of newborn cells (Figure 1B, Video 1). This was consistent with previous results in bulk cultures (Ondracka et al., 2018), confirming that our experimental setup for long-term live imaging faithfully reproduces culture growth. However, by measuring the diameter of the coenocytes, we observed that newborn cell release occurred with somewhat variable timing (Figure 1—figure supplement 1A, Video 1). Video 1 Download asset This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing. Download as MPEG-4 Download as WebM Download as Ogg Time lapse of synchronized cells of S. arctica obtained with epifluorescent microscopy. Time interval between frames is 20 min. The movie is played at 7fps. Four distinct cells can be seen undergoing a full life-cycle with the release of new born cells. Bar, 50 µm. Time-lapse imaging revealed that prior to cell release, the coenocytes darken along the periphery, and the dark front begins to advance towards the center (Video 2). Afterwards, we observed an abrupt internal morphological change in the coenocyte, when the front disappears. We termed this event ‘flip’ (Video 2). The flip occurred in all the coenocytes and can be reliably detected on kymographs (Figure 1C). Aligning individual coenocyte size traces to this specific common temporal marker, we observe that coenocytes stop growing approximately 9 hr before the flip (Figure 1D and E and Figure 1—figure supplement 1B). This shows that the growth stage and cellularization are temporally separated. The cellularization can be divided into a temporally less variable pre-flip phase (~9 hr) and a variable post-flip phase (Figure 1E and Figure 1—figure supplement 1B and C). Taken together, these results show that growth and cellularization form temporally distinct stages of the coenocytic cycle of S. arctica. This provides a temporal framework in which to characterize phenotypically distinct stages of cellularization. Video 2 Download asset This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing. Download as MPEG-4 Download as WebM Download as Ogg Time lapse of synchronized cells of S. arctica obtained with epifluorescent microscopy. Time interval between frames is 15 min. The movie is played at 7fps. Three cells can be seen undergoing flip prior to cell release whereas one cell undergo ‘flip’ without cell release during the course of the experiment. The exact timing of cell release is pointed out for each cell. Bar, 50 µm. The cortical actin network establishes sites of membrane invagination and generates a polarized layer of cells during cellularization To assess whether cellularization in S. arctica involves encapsulation of nuclei by plasma membranes, we imaged the plasma membrane using live time-lapse imaging in the presence of the lipid dye FM4-64 (Betz et al., 1996). We observed a rapid increase in FM4-64 intensity at the periphery of the coenocyte 30 min prior to flip (Figure 2A, panel II, Videos 3 and 4). The plasma membrane invagination sites formed at the periphery and progressed synchronously from the outside toward the center of the coenocyte, forming polarized, polyhedral cells (Figure 2A, panels II-V, Videos 3 and 4). Similar to cellularization in Drosophila coenocytes (Lecuit and Wieschaus, 2000; Lecuit et al., 2002), we observe different rates of plasma membrane invaginations throughout the process. During the first 10–15 min, the invaginations advance at a rate of about 0,3 µm/min (Figure 2—figure supplement 1A). This period is followed by a rapid phase where membrane invaginates at a rate of about 0,9 µm/min leading to flip (Figure 2—figure supplement 1A). Finally, following flip, cells lost their polyhedral shape and became round, suggesting that they were no longer attached to each other (Figure 2A, panel VI, Videos 3 and 4). Figure 2 with 1 supplement see all Download asset Open asset Actin cytoskeleton and plasma membrane dynamics during cellularization of S. arctica. (A) Dynamics of plasma membrane invaginations during cellularization. Live-cells, pre-grown for 58 hr, were stained with FM4-64 (10 µM) and imaged using epifluorescent microscopy with a 5 min interval. Bar, 50 µm. (B) Spatio-temporal organization of the actin cytoskeleton, nuclei and cells during cellularization. Synchronized cells of S. arctica, pre-grown for 48 hr, were fixed every hour for 14 hr and stained with phalloidin and DAPI to reveal cytoskeletal dynamics during cellularization. All cells were imaged using confocal microscopy. In panel IV, arrows indicate higher actin signal intensity on the internal side and that nuclei are localized close to the cortex indicating that the layer of cells is polarized. Bar, 50 µm. (C) Actin network is established prior to plasma membrane invaginations. Synchronized cells of S. arctica, pre-grown for 54 hr, were fixed every 2 hr for 10 hr and stained with both the membrane dye FM4-64FX and phalloidin. Arrows show sites of colocalization between both markers at the onset of plasma membrane invaginations. Bar, 50 µm. Video 3 Download asset This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing. Download as MPEG-4 Download as WebM Download as Ogg Time lapse of synchronized cells of S. arctica stained with the plasma membrane dye FM4-64 and obtained with epifluorescent microscopy. Time interval between frames is 5 min. The movie is played at 7fps. Plasma membrane invaginations can be seen occurring from the outside inwards for approximately 40 min prior to flip. Bar, 50 µm. Video 4 Download asset This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing. Download as MPEG-4 Download as WebM Download as Ogg Time lapse of plasma membrane dynamics during cellularization of S. arctica. In animal coenocytes, plasma membrane invagination is associated with dynamic organization of the actomyosin cytoskeleton (Mazumdar and Mazumdar, 2002). To investigate actin dynamics during cellularization, we took advantage of the timeline described above and imaged coenocytes that were fixed and stained for actin and nuclei (using phalloidin and DAPI, respectively) at different time points during cellularization (Figure 2B, Figure 2—figure supplement 1B and C). During growth, actin localized exclusively as small patches at the surface of the coenocyte (Figure 2—figure supplement 1B, panel I). Only at the onset of cellularization, multiple actin patches increased in size to form actin nodes (Figure 2B, panel I and II, and Figure 2—figure supplement 1B, panel II). This phase was followed by cortical compartmentalization surrounding the nuclei through gradual formation of an actin filament network solely at the cortex of the coenocyte (Figure 2B, panel III, and Figure 2—figure supplement 1B, panels III and IV). Following this cortical compartmentalization, a layer of cells was transiently formed by inward growth of the actin filaments from the cortex (Figure 2B, panels IV and Figure 2—figure supplement 1B, panels V and VI). During this stage, the actin signal intensity was uneven and higher on the internal side (Video 5), and nuclei were localized close to the cortex, indicating that cells are polarized (Figure 2B, panel IV, Video 5). These polarized cells progressively grew towards the center of the coenocyte to fill the cavity (Figure 2—figure supplement 1B, panel VII). After flip, similar to the plasma membrane organization mentioned above, the layer of cells was reorganized to form spherical cells (Figure 2B, panel VI, and Figure 2—figure supplement 1B, panel VIII). Video 5 Download asset This video cannot be played in place because your browser does support HTML5 video. You may still download the video for offline viewing. Download as MPEG-4 Download as WebM Download as Ogg Z-projection of the spatial organization of the polarized layer of cells during cellularization of S. arctica obtained with confocal microscopy. The cell is stained with phalloidin and DAPI. The movie is played at 7fps. A polarized layer of cells can be seen. Bar, 10 µm. To determine the order of actin organization and plasma membrane invaginations, we assessed the localization of both actin and plasma membrane in fixed samples using phalloidin and a fixable analog of FM4-64 (FM4-64FX). We found that the cortical actin network formed prior to the appearance of the membrane dye (Figure 2C, panel II). The intensity of FM4-64FX labeling also increased and co-localized with the underlying actin network around the time of plasma membrane invagination (Figure 2C, panels III and IV). This suggests that the cortical actin network determines the site of plasma membrane invagination. Finally, to determine the timing of cell wall formation, we stained the cells with calcofluor. We observed that labeling co-localized with the membrane marker FM4-64FX around individual cells in the post-flip stage (Figure 2—figure supplement 1D). This shows that the newborn cells already build the cell wall before the release, as was suggested previously in other Sphaeroforma species (Marshall and Berbee, 2013). In early insect embryos, cellularization of the syncytial blastoderm occurs through actin-dependent invagination of the plasma membrane. Here, we demonstrate that the cellularization of the ichthyosporean coenocyte also involves active actin reorganization and membrane invagination (slow and fast phases of membrane invagination rates) at the site of actin cytoskeleton. Additionally, this results in the transient formation of a polarized layer of cells with an internal cavity that morphologically resembles simple epithelial structures. Cellularization is associated with extensive sequential transcriptional waves and is associated with evolutionarily younger transcripts To gain insight into the regulation of the cellularization of S. arctica, we sequenced the expressed mRNAs of synchronized cultures with high temporal resolution, and comprehensively analyzed the dynamics of transcription, alternative splicing, and long intergenic non-coding RNAs (lincRNAs). Because the published genome assembly of S. arctica (Grau-Bové et al., 2017) was fragmented and likely resulted in incomplete gene models, we first re-sequenced the genome combining the Illumina technology with the PacBio long read sequencing technology. The final assembly sequences comprised 142,721,209 bp, and the metrics were greatly improved compared to the previous assembly (Grau-Bové et al., 2017) (Figure 3—source data 1). Ab initio gene annotation resulted in the discovery of novel ORFs due to the absence of repetitive regions in the previous assembly. RNA-seq data was used to improve the ORF prediction, to define the 5’ and 3’ untranslated regions, and to discover lincRNAs. In total, 33,682 protein coding genes and 1071 lincRNAs were predicted using this combined pipeline (see Materilas and methods). This final transcriptome assembly was used as the reference transcriptome for further analysis. To perform the time-resolved transcriptomics, we isolated and sequenced mRNA from two independent synchronized cultures at 6 hr time intervals during the entire coenocytic cycle, encompassing time points from early 4-nuclei stage throughout growth and cellularization stages until the release of newborn cells (Figure 3—figure supplement 1A). We first analyzed transcript abundance during the time series. The majority of the transcriptome (20,196 out of 34,753 predicted protein-coding and lincRNA genes) was transcribed at very low levels (mean expression throughout the time course <0.5 tpm [transcripts per million]) and were removed for the subsequent clustering analysis. Clustering of transcript abundance data from both biological replicates revealed a clear separation between the transcriptomes of the growth stage (12 hr to 42 hr time points) and the cellularization stage (48 hr to 66 hr time points) (Figure 3A). Furthermore, we observed that the transcriptome patterns between replicates 1 and 2 were shifted by 6 hr from 48 hr onwards (Figure 3A), presumably due to differences in temperature and conditions influencing the kinetics of the coenocytic cycle. We thus adjusted the time of the second replicate by 6 hr according to the clustering results, although we emphasize that the clustering analysis did not depend on time. Among the expressed transcripts (defined as mean expression levels higher than 0.5 tpm across all samples; in total 13,542 coding genes and 1015 lincRNAs), consensus clustering using Clust (Abu-Jamous and Kelly, 2018) extracted 9 clusters of co-expressed transcripts with a total of 4441 protein coding genes (Figure 3B), while the rest of the transcripts were not assigned to any co-expression cluster. The assigned cluster membership was robust to using either of the replicates or averaging (Figure 3—figure supplement 1B). Visualization by heatmap and t-SNE plot separated the clusters into two meta-clusters containing the growth stage (clusters 1–3, totaling 2197 genes) and cellularization stage (clusters 4–9, totaling 2314 genes) clusters (Figure 3C and D). Among the cellularization clusters, we obtain three main clusters of genes sharply activated during cellularization (clusters 5, 6 and 7), which contain transcripts expressed at early, mid and late cellularization, respectively (Figure 3B and D). Altogether, this shows that cellularization is associated with extensive sharp transcriptional activation in multiple temporal waves, in total affecting 17% of the expressed transcriptome. Figure 3 with 1 supplement see all Download asset Open asset Transcriptional dynamics across the S. arctica life cycle. (A) Hierarchical clustering of time point samples by Euclidian distance of spearman correlation coefficient. The sample T48_2 is missing due to technical reason. (B) A heatmap of 4441 coding genes that were clustered into nine clusters. (C) A t-SNE plot of clustered genes. (D) Mean expression profile of each gene expression cluster. Figure 3—source data 1 Metrics of the Sphaeroforma arctica genome assemblies. https://cdn.elifesciences.org/articles/49801/elife-49801-fig3-data1-v4.pdf Download elife-49801-fig3-data1-v4.pdf Figure 3—source data 2 Transcript abundance (in tpm) of all the S. arctica transcripts. Used in Figures 3, 4 and 5. https://cdn.elifesciences.org/articles/49801/elife-49801-fig3-data2-v4.txt Download elife-49801-fig3-data2-v4.txt Figure 3—source data 3 Normalized transcript abundance of the expressed S. arctica transcripts (mean tpm > 0.5). Used in Figures 3, 4 and 5. https://cdn.elifesciences.org/articles/49801/elife-49801-fig3-data3-v4.txt Download elife-49801-fig3-data3-v4.txt Figure 3—source data 4 Table of transcripts per cluster membership. https://cdn.elifesciences.org/articles/49801/elife-49801-fig3-data4-v4.txt Download elife-49801-fig3-data4-v4.txt In parallel to transcript abundance dynamics, we also assessed alternative splicing (AS) across the time series. This analysis identified 2022 genes affected by intron retention (12.9% of all intron-bearing genes, totaling 4310 introns), 914 by exon skipping (12.3% of genes with >2 exons, 1206 exons) and 44 with mutually exclusive exons (0.7% of genes with >3 exons, involving 118 exon pairs) in all samples (Figure 3—figure supplement 1C and D). Overall, neither the number of AS events nor the number of genes affected vary dramatically along the S. arctica growth cycle (Figure 3—figure supplement 1E-J). Analysis of AS events over time did not yield any discernible global dynamics, although we found a small number of events differentially present between the growth and cellularization stages (Figure 3—figure supplement 1K). Interestingly, skipped exons were more likely to be in-frame (38.63%, compared to 30.33% of in-frame exons in genes with >2 exons, p=4.34e-05, Fisher’s exact test) and yield non-truncated transcripts, a phenomenon commonly observed in animal transcriptomes but not in transcriptomes of other unicellular eukaryotes (Grau-Bové et al., 2018). The transcripts affected by such in-frame exon skipping events are enriched in biological processes such as regulation of multicellular organismal processes, assembly of the focal adhesion complex and cell growth (Figure 3—figure supplement 1L). Overall, although pervasive, alternative splicing likely does not play a major role in regulation of the coenocytic cycle and cellularization of S. arctica. Additionally, we analyzed the dynamics of lincRNA expression. Overall, lincRNAs represent ~3% of total transcript abundance across the time series (Figure 4—figure supplement 1A). Among the long non-coding RNAs, 70 lincRNA transcripts clustered with coding genes into temporally co-expressed clusters (Figure 4A). Sequence homology searches revealed that 24 of the S. arctica lncRNAs were conserved in distantly related ichthyosporean species such as Creolimax fragrantissima, Pirum gemmata and Abeoforma whisleri (estimated to have diverged ~500 million years ago; Parfrey et al., 2011). This is a remarkable depth of conservation, since animal lincRNAs are not conserved between animal phyla (Bråte et al., 2015; Gaiti et al., 2015; Hezroni et al., 2015). Other lincRNAs were either specific to S. arctica (511) or conserved only in closely related Sphaeroforma species (536). Comparison of lincRNAs by degree of conservation showed no notable difference in GC content (Figure 4—figure supplement 1B) or dynamics of expression during the coenocytic cycle (Figure 4—figure supplement 1C). However, we found that conserved lincRNAs were on average longer than non-conserved ones (Figure 4—figure supplement 1D) and, strikingly, expressed at much higher levels (Figure 4B). Among the 24 deeply conserved lincRNAs, five clustered in the temporally co-expressed clusters, which is higher than expected by chance (p=0.0166, Fisher’s exact test). Among these, three belonged to the cellularization clusters, including, lincRNA asmbl_31839, which has a remarkably high sequence similarity with its homologs from other ichthyosporeans (Figure 4—figure supplement 1E). Furthermore, its transcriptional regulation is independent of the transcriptional regulation of its neighboring coding genes (located within 3 kb) (Figure 4C). In summary, we discovered deeply conserved lincRNAs that are expressed at high levels and are transcriptionally activated during cellularization. Figure 4 with 2 supplements see all Download asset Open asset Dynamics of lincRNAs, alternative splicing, and gene phylostrata. (A) A heatmap of 70 long intergenic non-coding RNAs (lincRNAs) that co-cluster with coding genes. (B) mean expression level of the lincRNAs, binned by degree of conservation. (C) Expression of the conserved lincRNA, lincRNA: asmbl_31839 and the two coding genes located immediately upstream and downstream of it. (D) A phylogenetic tree indicating the 10 defined phylostrata, and a pie chart of fractions of all expressed genes per phylostratum. (E) Coefficient of variance of gene expression across the S. arctica life cycle, binned by phylostratum. (F) Transcriptome age index (TAI) of gene expression for the S. arctica life cycle. Higher TAI represents a younger transcriptome. Figure 4—source data 1 Table of Blast e-values of orthologs of S. arctica lincRNAs in other ichthyosporean species. https://cdn.elifesciences.org/articles/49801/elife-49801-fig4-data1-v4.txt Download elife-49801-fig4-data1-v4.txt Figure 4—source data 2 List of transcriptomes used for lincRNA conservation analysis. https://cdn.elifesciences.org/articles/49801/elife-49801-fig4-data2-v4.pdf Download elife-49801-fig4-data2-v4.pdf Figure 4—source data 3 Table of gene orthogroups. Used in Figures 4 and 5. https://cdn.elifesciences.org/articles/49801/elife-49801-fig4-data3-v4.csv Download elife-49801-fig4-data3-v4.csv Figure 4—source data 4 List of transcriptomes used for generation of orthogroups and phylostratigraphy analysis. https://cdn.elifesciences.org/articles/49801/elife-49801-fig4-data4-v4.pdf Download elife-49801-fig4-data4-v4.pdf Figure 4—source data 5 Table of S. arctica transcripts per gene stratum and orthogroup membership. Gene ages range from 1 (S. arctica-specific) to 10 (pan-eukaryotic). Only transcripts that were classified in an orthogroup are listed; the remainder of the transcripts were S. arctica-specific. https://cdn.elifesciences.org/articles/49801/elife-49801-fig4-data5-v4.csv Download elife-49801-fig4-data5-v4.csv Finally, to assess the evolutionary origins of the co-expression clusters, we used a phylostratigraphic analysis to classify genes into evolutionary gene age groups. We carried out orthology analysis of the predicted S. arctica proteome along with 30 representative species from the eukaryotic tree of life to identify ‘orthogroups’ (i.e. groups of putative orthologs between species). S. arctica protein-coding genes clustered in 6149 orthogroups representing 12,527 genes; the rest of the genes did not have an ortholog outside S. arctica. Next, we inferred the age of each gene using Dollo parsimony (Csurös, 2010) to classify them into phylostrata (sets of genes with the same phylogenetic origin) (Figure 4D). Analysis of gene expression by phylostrata revealed a trend toward more variable expression throughout the coenocytic cycle in younger genes (Figure 4E), although their mean expression levels were lower (Figure 4—figure supplement 2A). Such a correlation has been observed in animal development, where developmentally regulated genes tend to be of younger origin (Domazet-Lošo and Tautz, 2010). Analysis of enrichment of gene phylostrata in each gene expression cluster (Figure 4—figure supplement 2B) showed that the growth clusters are enriched for pan-eukaryotic genes. In contrast, we find that the cellularization clusters were enriched for younger genes (Figure 4—figure supplement 2B and C). Importantly, we found that genes with ichthyosporean origins were significantly enriched in all three of the largest cellularization clusters (Figure 4—figure supplement 2B and C). Finally, computing the transcriptome age index (Domazet-Lošo and Tautz, 2010) to assess the overall age of the transcriptome across the life cycle revealed a statistically significant hourglass pattern (p=4.4 * 10−5 and p=0.02 for replicates 1 and 2, Reductive Hourglass test; Drost et al., 2015), with older genes expressed at later stages of growth, and younger genes expressed during early growth and cellu