Spatial profiling of early primate gastrulation in utero.

Gastrulation controls the emergence of cellular diversity and axis patterning in the early embryo. In mammals, this transformation is orchestrated by dynamic signalling centres at the interface of embryonic and extraembryonic tissues1-3. Elucidating the molecular framework of axis formation in vivo is fundamental for our understanding of human development4-6 and to advance stem-cell-based regenerative approaches7. Here we illuminate early gastrulation of marmoset embryos in utero using spatial transcriptomics and stem-cell-based embryo models. Gaussian process regression-based 3D transcriptomes delineate the emergence of the anterior visceral endoderm, which is hallmarked by conserved (HHEX, LEFTY2, LHX1) and primate-specific (POSTN, SDC4, FZD5) factors. WNT signalling spatially coordinates the formation of the primitive streak in the embryonic disc and is counteracted by SFRP1 and SFRP2 to sustain pluripotency in the anterior domain. Amnion specification occurs at the boundaries of the embryonic disc through ID1, ID2 and ID3 in response to BMP signalling, providing a developmental rationale for amnion differentiation of primate pluripotent stem cells (PSCs). Spatial identity mapping demonstrates that primed marmoset PSCs exhibit the highest similarity to the anterior embryonic disc, whereas naive PSCs resemble the preimplantation epiblast. Our 3D transcriptome models reveal the molecular code of lineage specification in the primate embryo and provide an in vivo reference to decipher human development.

Microgel culture and spatial identity mapping elucidate the signalling requirements for primate epiblast and amnion formation.

The early specification and rapid growth of extraembryonic membranes are distinctive hallmarks of primate embryogenesis. These complex tasks are resolved through an intricate combination of signals controlling the induction of extraembryonic lineages and, at the same time, safeguarding the pluripotent epiblast. Here, we delineate the signals orchestrating primate epiblast and amnion identity. We encapsulated marmoset pluripotent stem cells into agarose microgels and identified culture conditions for the development of epiblast- and amnion-spheroids. Spatial identity mapping authenticated spheroids generated in vitro by comparison with marmoset embryos in vivo. We leveraged the microgel system to functionally interrogate the signalling environment of the post-implantation primate embryo. Single-cell profiling of the resulting spheroids demonstrated that activin/nodal signalling is required for embryonic lineage identity. BMP4 promoted amnion formation and maturation, which was counteracted by FGF signalling. Our combination of microgel culture, single-cell profiling and spatial identity mapping provides a powerful approach to decipher the essential cues for embryonic and extraembryonic lineage formation in primate embryogenesis.

Time-series transcriptomics reveals a BBX32-directed control of acclimation to high light in mature Arabidopsis leaves.

The photosynthetic capacity of mature leaves increases after several days' exposure to constant or intermittent episodes of high light (HL) and is manifested primarily as changes in chloroplast physiology. How this chloroplast-level acclimation to HL is initiated and controlled is unknown. From expanded Arabidopsis leaves, we determined HL-dependent changes in transcript abundance of 3844 genes in a 0-6 h time-series transcriptomics experiment. It was hypothesized that among such genes were those that contribute to the initiation of HL acclimation. By focusing on differentially expressed transcription (co-)factor genes and applying dynamic statistical modelling to the temporal transcriptomics data, a regulatory network of 47 predominantly photoreceptor-regulated transcription (co-)factor genes was inferred. The most connected gene in this network was B-BOX DOMAIN CONTAINING PROTEIN32 (BBX32). Plants overexpressing BBX32 were strongly impaired in acclimation to HL and displayed perturbed expression of photosynthesis-associated genes under LL and after exposure to HL. These observations led to demonstrating that as well as regulation of chloroplast-level acclimation by BBX32, CRYPTOCHROME1, LONG HYPOCOTYL5, CONSTITUTIVELY PHOTOMORPHOGENIC1 and SUPPRESSOR OF PHYA-105 are important. In addition, the BBX32-centric gene regulatory network provides a view of the transcriptional control of acclimation in mature leaves distinct from other photoreceptor-regulated processes, such as seedling photomorphogenesis.

Branch-recombinant Gaussian processes for analysis of perturbations in biological time series.

Motivation: A common class of behaviour encountered in the biological sciences involves branching and recombination. During branching, a statistical process bifurcates resulting in two or more potentially correlated processes that may undergo further branching; the contrary is true during recombination, where two or more statistical processes converge. A key objective is to identify the time of this bifurcation (branch or recombination time) from time series measurements, e.g. by comparing a control time series with perturbed time series. Gaussian processes (GPs) represent an ideal framework for such analysis, allowing for nonlinear regression that includes a rigorous treatment of uncertainty. Currently, however, GP models only exist for two-branch systems. Here, we highlight how arbitrarily complex branching processes can be built using the correct composition of covariance functions within a GP framework, thus outlining a general framework for the treatment of branching and recombination in the form of branch-recombinant Gaussian processes (B-RGPs). Results: We first benchmark the performance of B-RGPs compared to a variety of existing regression approaches, and demonstrate robustness to model misspecification. B-RGPs are then used to investigate the branching patterns of Arabidopsis thaliana gene expression following inoculation with the hemibotrophic bacteria, Pseudomonas syringae DC3000, and a disarmed mutant strain, hrpA. By grouping genes according to the number of branches, we could naturally separate out genes involved in basal immune response from those subverted by the virulent strain, and show enrichment for targets of pathogen protein effectors. Finally, we identify two early branching genes WRKY11 and WRKY17, and show that genes that branched at similar times to WRKY11/17 were enriched for W-box binding motifs, and overrepresented for genes differentially expressed in WRKY11/17 knockouts, suggesting that branch time could be used for identifying direct and indirect binding targets of key transcription factors. Availability and implementation: https://github.com/cap76/BranchingGPs. Supplementary information: Supplementary data are available at Bioinformatics online.

Time-Series Transcriptomics Reveals That AGAMOUS-LIKE22 Affects Primary Metabolism and Developmental Processes in Drought-Stressed Arabidopsis.

In Arabidopsis thaliana, changes in metabolism and gene expression drive increased drought tolerance and initiate diverse drought avoidance and escape responses. To address regulatory processes that link these responses, we set out to identify genes that govern early responses to drought. To do this, a high-resolution time series transcriptomics data set was produced, coupled with detailed physiological and metabolic analyses of plants subjected to a slow transition from well-watered to drought conditions. A total of 1815 drought-responsive differentially expressed genes were identified. The early changes in gene expression coincided with a drop in carbon assimilation, and only in the late stages with an increase in foliar abscisic acid content. To identify gene regulatory networks (GRNs) mediating the transition between the early and late stages of drought, we used Bayesian network modeling of differentially expressed transcription factor (TF) genes. This approach identified AGAMOUS-LIKE22 (AGL22), as key hub gene in a TF GRN. It has previously been shown that AGL22 is involved in the transition from vegetative state to flowering but here we show that AGL22 expression influences steady state photosynthetic rates and lifetime water use. This suggests that AGL22 uniquely regulates a transcriptional network during drought stress, linking changes in primary metabolism and the initiation of stress responses.

Comparative transcriptome analysis in Arabidopsis ein2/ore3 and ahk3/ore12 mutants during dark-induced leaf senescence.

Leaf senescence involves degenerative but active biological processes that require balanced regulation of pro- and anti-senescing activities. Ethylene and cytokinin are major antagonistic regulatory hormones that control the timing and progression rate of leaf senescence. To identify the roles of these hormones in the regulation of leaf senescence in Arabidopsis, global gene expression profiles in detached leaves of the wild type, an ethylene-insensitive mutant (ein2/ore3), and a constitutive cytokinin response mutant (ahk3/ore12) were investigated during dark-induced leaf senescence. Comparative transcriptome analyses revealed that genes involved in oxidative or salt stress response were preferentially altered in the ein2/ore3 mutant, whereas genes involved in ribosome biogenesis were affected in the ahk3/ore12 mutant during dark-induced leaf senescence. Similar results were also obtained for developmental senescence. Through extensive molecular and physiological analyses in ein2/ore3 and ahk3/ore12 during dark-induced leaf senescence, together with responses when treated with cytokinin and ethylene inhibitor, we conclude that ethylene acts as a senescence-promoting factor via the transcriptional regulation of stress-related responses, whereas cytokinin acts as an anti-senescing agent by maintaining cellular activities and preserving the translational machinery. These findings provide new insights into how plants utilize two antagonistic hormones, ethylene and cytokinin, to regulate the molecular programming of leaf senescence.

Transcriptional Dynamics Driving MAMP-Triggered Immunity and Pathogen Effector-Mediated Immunosuppression in Arabidopsis Leaves Following Infection with Pseudomonas syringae pv tomato DC3000.

Transcriptional reprogramming is integral to effective plant defense. Pathogen effectors act transcriptionally and posttranscriptionally to suppress defense responses. A major challenge to understanding disease and defense responses is discriminating between transcriptional reprogramming associated with microbial-associated molecular pattern (MAMP)-triggered immunity (MTI) and that orchestrated by effectors. A high-resolution time course of genome-wide expression changes following challenge with Pseudomonas syringae pv tomato DC3000 and the nonpathogenic mutant strain DC3000hrpA- allowed us to establish causal links between the activities of pathogen effectors and suppression of MTI and infer with high confidence a range of processes specifically targeted by effectors. Analysis of this information-rich data set with a range of computational tools provided insights into the earliest transcriptional events triggered by effector delivery, regulatory mechanisms recruited, and biological processes targeted. We show that the majority of genes contributing to disease or defense are induced within 6 h postinfection, significantly before pathogen multiplication. Suppression of chloroplast-associated genes is a rapid MAMP-triggered defense response, and suppression of genes involved in chromatin assembly and induction of ubiquitin-related genes coincide with pathogen-induced abscisic acid accumulation. Specific combinations of promoter motifs are engaged in fine-tuning the MTI response and active transcriptional suppression at specific promoter configurations by P. syringae.

Inferring the perturbation time from biological time course data.

MOTIVATION: Time course data are often used to study the changes to a biological process after perturbation. Statistical methods have been developed to determine whether such a perturbation induces changes over time, e.g. comparing a perturbed and unperturbed time course dataset to uncover differences. However, existing methods do not provide a principled statistical approach to identify the specific time when the two time course datasets first begin to diverge after a perturbation; we call this the perturbation time. Estimation of the perturbation time for different variables in a biological process allows us to identify the sequence of events following a perturbation and therefore provides valuable insights into likely causal relationships. RESULTS: We propose a Bayesian method to infer the perturbation time given time course data from a wild-type and perturbed system. We use a non-parametric approach based on Gaussian Process regression. We derive a probabilistic model of noise-corrupted and replicated time course data coming from the same profile before the perturbation time and diverging after the perturbation time. The likelihood function can be worked out exactly for this model and the posterior distribution of the perturbation time is obtained by a simple histogram approach, without recourse to complex approximate inference algorithms. We validate the method on simulated data and apply it to study the transcriptional change occurring in Arabidopsis following inoculation with Pseudomonas syringae pv. tomato DC3000 versus the disarmed strain DC3000hrpA AVAILABILITY AND IMPLEMENTATION: : An R package, DEtime, implementing the method is available at https://github.com/ManchesterBioinference/DEtime along with the data and code required to reproduce all the results. CONTACT: Jing.Yang@manchester.ac.uk or Magnus.Rattray@manchester.ac.uk SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Inferring orthologous gene regulatory networks using interspecies data fusion.

MOTIVATION: The ability to jointly learn gene regulatory networks (GRNs) in, or leverage GRNs between related species would allow the vast amount of legacy data obtained in model organisms to inform the GRNs of more complex, or economically or medically relevant counterparts. Examples include transferring information from Arabidopsis thaliana into related crop species for food security purposes, or from mice into humans for medical applications. Here we develop two related Bayesian approaches to network inference that allow GRNs to be jointly inferred in, or leveraged between, several related species: in one framework, network information is directly propagated between species; in the second hierarchical approach, network information is propagated via an unobserved 'hypernetwork'. In both frameworks, information about network similarity is captured via graph kernels, with the networks additionally informed by species-specific time series gene expression data, when available, using Gaussian processes to model the dynamics of gene expression. RESULTS: Results on in silico benchmarks demonstrate that joint inference, and leveraging of known networks between species, offers better accuracy than standalone inference. The direct propagation of network information via the non-hierarchical framework is more appropriate when there are relatively few species, while the hierarchical approach is better suited when there are many species. Both methods are robust to small amounts of mislabelling of orthologues. Finally, the use of Saccharomyces cerevisiae data and networks to inform inference of networks in the budding yeast Schizosaccharomyces pombe predicts a novel role in cell cycle regulation for Gas1 (SPAC19B12.02c), a 1,3-beta-glucanosyltransferase. AVAILABILITY AND IMPLEMENTATION: MATLAB code is available from http://go.warwick.ac.uk/systemsbiology/software/.

CSI: a nonparametric Bayesian approach to network inference from multiple perturbed time series gene expression data.

Here we introduce the causal structure identification (CSI) package, a Gaussian process based approach to inferring gene regulatory networks (GRNs) from multiple time series data. The standard CSI approach infers a single GRN via joint learning from multiple time series datasets; the hierarchical approach (HCSI) infers a separate GRN for each dataset, albeit with the networks constrained to favor similar structures, allowing for the identification of context specific networks. The software is implemented in MATLAB and includes a graphical user interface (GUI) for user friendly inference. Finally the GUI can be connected to high performance computer clusters to facilitate analysis of large genomic datasets.

Correcting for link loss in causal network inference caused by regulator interference.

MOTIVATION: There are a number of algorithms to infer causal regulatory networks from time series (gene expression) data. Here we analyse the phenomena of regulator interference, where regulators with similar dynamics mutually suppress both the probability of regulating a target and the associated link strength; for instance, interference between two identical strong regulators reduces link probabilities by ∼50%. RESULTS: We construct a robust method to define an interference-corrected causal network based on an analysis of the conditional link probabilities that recovers links lost through interference. On a large real network (Streptomyces coelicolor, phosphate depletion), we demonstrate that significant interference can occur between regulators with a correlation as low as 0.865, losing an estimated 34% of links by interference. However, levels of interference cannot be predicted from the correlation between regulators alone and are data specific. Validating against known networks, we show that high numbers of functional links are lost by regulator interference. Performance against other methods on DREAM4 data is excellent. AVAILABILITY AND IMPLEMENTATION: The method is implemented in R and is publicly available as the NIACS package at http://www2.warwick.ac.uk/fac/sci/systemsbiology/research/software.

Arabidopsis defense against Botrytis cinerea: chronology and regulation deciphered by high-resolution temporal transcriptomic analysis.

Transcriptional reprogramming forms a major part of a plant's response to pathogen infection. Many individual components and pathways operating during plant defense have been identified, but our knowledge of how these different components interact is still rudimentary. We generated a high-resolution time series of gene expression profiles from a single Arabidopsis thaliana leaf during infection by the necrotrophic fungal pathogen Botrytis cinerea. Approximately one-third of the Arabidopsis genome is differentially expressed during the first 48 h after infection, with the majority of changes in gene expression occurring before significant lesion development. We used computational tools to obtain a detailed chronology of the defense response against B. cinerea, highlighting the times at which signaling and metabolic processes change, and identify transcription factor families operating at different times after infection. Motif enrichment and network inference predicted regulatory interactions, and testing of one such prediction identified a role for TGA3 in defense against necrotrophic pathogens. These data provide an unprecedented level of detail about transcriptional changes during a defense response and are suited to systems biology analyses to generate predictive models of the gene regulatory networks mediating the Arabidopsis response to B. cinerea.

A local regulatory network around three NAC transcription factors in stress responses and senescence in Arabidopsis leaves.

A model is presented describing the gene regulatory network surrounding three similar NAC transcription factors that have roles in Arabidopsis leaf senescence and stress responses. ANAC019, ANAC055 and ANAC072 belong to the same clade of NAC domain genes and have overlapping expression patterns. A combination of promoter DNA/protein interactions identified using yeast 1-hybrid analysis and modelling using gene expression time course data has been applied to predict the regulatory network upstream of these genes. Similarities and divergence in regulation during a variety of stress responses are predicted by different combinations of upstream transcription factors binding and also by the modelling. Mutant analysis with potential upstream genes was used to test and confirm some of the predicted interactions. Gene expression analysis in mutants of ANAC019 and ANAC055 at different times during leaf senescence has revealed a distinctly different role for each of these genes. Yeast 1-hybrid analysis is shown to be a valuable tool that can distinguish clades of binding proteins and be used to test and quantify protein binding to predicted promoter motifs.

Modelling transcriptional networks in leaf senescence.

The process of leaf senescence is induced by an extensive range of developmental and environmental signals and controlled by multiple, cross-linking pathways, many of which overlap with plant stress-response signals. Elucidation of this complex regulation requires a step beyond a traditional one-gene-at-a-time analysis. Application of a more global analysis using statistical and mathematical tools of systems biology is an approach that is being applied to address this problem. A variety of modelling methods applicable to the analysis of current and future senescence data are reviewed and discussed using some senescence-specific examples. Network modelling with a senescence transcriptome time course followed by testing predictions with gene-expression data illustrates the application of systems biology tools.

High-resolution temporal profiling of transcripts during Arabidopsis leaf senescence reveals a distinct chronology of processes and regulation.

Leaf senescence is an essential developmental process that impacts dramatically on crop yields and involves altered regulation of thousands of genes and many metabolic and signaling pathways, resulting in major changes in the leaf. The regulation of senescence is complex, and although senescence regulatory genes have been characterized, there is little information on how these function in the global control of the process. We used microarray analysis to obtain a high-resolution time-course profile of gene expression during development of a single leaf over a 3-week period to senescence. A complex experimental design approach and a combination of methods were used to extract high-quality replicated data and to identify differentially expressed genes. The multiple time points enable the use of highly informative clustering to reveal distinct time points at which signaling and metabolic pathways change. Analysis of motif enrichment, as well as comparison of transcription factor (TF) families showing altered expression over the time course, identify clear groups of TFs active at different stages of leaf development and senescence. These data enable connection of metabolic processes, signaling pathways, and specific TF activity, which will underpin the development of network models to elucidate the process of senescence.

Age-related and species-specific methylation changes in the protein-coding marmoset sperm epigenome.

The sperm epigenome is thought to affect the developmental programming of the resulting embryo, influencing health and disease in later life. Age-related methylation changes in the sperm of old fathers may mediate the increased risks for reproductive and offspring medical problems. The impact of paternal age on sperm methylation has been extensively studied in humans and, to a lesser extent, in rodents and cattle. Here, we performed a comparative analysis of paternal age effects on protein-coding genes in the human and marmoset sperm methylomes. The marmoset has gained growing importance as a non-human primate model of aging and age-related diseases. Using reduced representation bisulfite sequencing, we identified age-related differentially methylated transcription start site (ageTSS) regions in 204 marmoset and 27 human genes. The direction of methylation changes was the opposite, increasing with age in marmosets and decreasing in humans. None of the identified ageTSS was differentially methylated in both species. Although the average methylation levels of all TSS regions were highly correlated between marmosets and humans, with the majority of TSS being hypomethylated in sperm, more than 300 protein-coding genes were endowed with species-specifically (hypo)methylated TSS. Several genes of the glycosphingolipid (GSL) biosynthesis pathway, which plays a role in embryonic stem cell differentiation and regulation of development, were hypomethylated (<5%) in human and fully methylated (>95%) in marmoset sperm. The expression levels and patterns of defined sets of GSL genes differed considerably between human and marmoset pre-implantation embryo stages and blastocyst tissues, respectively.

Nonparametric Bayesian inference for perturbed and orthologous gene regulatory networks.

MOTIVATION: The generation of time series transcriptomic datasets collected under multiple experimental conditions has proven to be a powerful approach for disentangling complex biological processes, allowing for the reverse engineering of gene regulatory networks (GRNs). Most methods for reverse engineering GRNs from multiple datasets assume that each of the time series were generated from networks with identical topology. In this study, we outline a hierarchical, non-parametric Bayesian approach for reverse engineering GRNs using multiple time series that can be applied in a number of novel situations including: (i) where different, but overlapping sets of transcription factors are expected to bind in the different experimental conditions; that is, where switching events could potentially arise under the different treatments and (ii) for inference in evolutionary related species in which orthologous GRNs exist. More generally, the method can be used to identify context-specific regulation by leveraging time series gene expression data alongside methods that can identify putative lists of transcription factors or transcription factor targets. RESULTS: The hierarchical inference outperforms related (but non-hierarchical) approaches when the networks used to generate the data were identical, and performs comparably even when the networks used to generate data were independent. The method was subsequently used alongside yeast one hybrid and microarray time series data to infer potential transcriptional switches in Arabidopsis thaliana response to stress. The results confirm previous biological studies and allow for additional insights into gene regulation under various abiotic stresses. AVAILABILITY: The methods outlined in this article have been implemented in Matlab and are available on request.

Modeling meiotic chromosomes indicates a size dependent contribution of telomere clustering and chromosome rigidity to homologue juxtaposition.

Meiosis is the cell division that halves the genetic component of diploid cells to form gametes or spores. To achieve this, meiotic cells undergo a radical spatial reorganisation of chromosomes. This reorganisation is a prerequisite for the pairing of parental homologous chromosomes and the reductional division, which halves the number of chromosomes in daughter cells. Of particular note is the change from a centromere clustered layout (Rabl configuration) to a telomere clustered conformation (bouquet stage). The contribution of the bouquet structure to homologous chromosome pairing is uncertain. We have developed a new in silico model to represent the chromosomes of Saccharomyces cerevisiae in space, based on a worm-like chain model constrained by attachment to the nuclear envelope and clustering forces. We have asked how these constraints could influence chromosome layout, with particular regard to the juxtaposition of homologous chromosomes and potential nonallelic, ectopic, interactions. The data support the view that the bouquet may be sufficient to bring short chromosomes together, but the contribution to long chromosomes is less. We also find that persistence length is critical to how much influence the bouquet structure could have, both on pairing of homologues and avoiding contacts with heterologues. This work represents an important development in computer modeling of chromosomes, and suggests new explanations for why elucidating the functional significance of the bouquet by genetics has been so difficult.

Rat post-implantation epiblast-derived pluripotent stem cells produce functional germ cells.

In mammals, pluripotent cells transit through a continuum of distinct molecular and functional states en route to initiating lineage specification. Capturing pluripotent stem cells (PSCs) mirroring in vivo pluripotent states provides accessible in vitro models to study the pluripotency program and mechanisms underlying lineage restriction. Here, we develop optimal culture conditions to derive and propagate post-implantation epiblast-derived PSCs (EpiSCs) in rats, a valuable model for biomedical research. We show that rat EpiSCs (rEpiSCs) can be reset toward the naive pluripotent state with exogenous Klf4, albeit not with the other five candidate genes (Nanog, Klf2, Esrrb, Tfcp2l1, and Tbx3) effective in mice. Finally, we demonstrate that rat EpiSCs retain competency to produce authentic primordial germ cell-like cells that undergo functional gametogenesis leading to the birth of viable offspring. Our findings in the rat model uncover principles underpinning pluripotency and germline competency across species.

Origin and segregation of the human germline.

Human germline-soma segregation occurs during weeks 2-3 in gastrulating embryos. Although direct studies are hindered, here, we investigate the dynamics of human primordial germ cell (PGCs) specification using in vitro models with temporally resolved single-cell transcriptomics and in-depth characterisation using in vivo datasets from human and nonhuman primates, including a 3D marmoset reference atlas. We elucidate the molecular signature for the transient gain of competence for germ cell fate during peri-implantation epiblast development. Furthermore, we show that both the PGCs and amnion arise from transcriptionally similar TFAP2A-positive progenitors at the posterior end of the embryo. Notably, genetic loss of function experiments shows that TFAP2A is crucial for initiating the PGC fate without detectably affecting the amnion and is subsequently replaced by TFAP2C as an essential component of the genetic network for PGC fate. Accordingly, amniotic cells continue to emerge from the progenitors in the posterior epiblast, but importantly, this is also a source of nascent PGCs.