Abstract Epithelial-to-mesenchymal transition (EMT) is activated to equip cells with the capacity to adapt to and escape hostile conditions. While EMT is required for cancer progression, its role in breast cancer initiation remains elusive. Given the basal-like phenotype of breast cancers arising in female carriers of germline BRCA1 pathogenic variants ( BRCA1 carriers), we hypothesized that enhanced EMT susceptibility underlies precancerous initiation in mammary epithelium. Perturbation of patient-derived normal mammary organoids from BRCA1 carriers and non-carriers with inflammatory cytokines induced copy number variations (CNV) and the acquisition of oncogenic mutations in both groups. However, in organoids derived from BRCA1 carriers, cytokine exposure induced morphological, transcriptomic, and functional EMT, accompanied by a transition to basal-like phenotype. Concomitant DNA damage accumulation in organoids from BRCA1 -carriers demonstrated PARP inhibitor sensitivity. EMT-primed states were identified in a subpopulation of normal mammary epithelium from BRCA1 carriers. We demonstrate the utility of patient-derived normal BRCA1 heterozygous mammary organoids to reveal a plastic, high-risk epithelial state that is associated with a transient, targetable vulnerability.

Summary The BMP and TGFβ signaling pathways control cellular fate decisions in diverse biological contexts, often playing opposing roles. Despite extensive knowledge of these pathways, understanding how cells respond to environments containing these opposing cues remains a challenge. Here, we systematically analyze the activation of these pathways under combinatorial signaling environments. We find that TGFβ ligands inhibit BMP signaling, while BMP ligands enhance TGFβ transcriptional response across concentrations, ligand variants, and cell types. This asymmetric crosstalk results in the activation of a TGFβ-biased transcriptional response, even under mixed signaling conditions, effectively reducing signal ambiguity, with implications for processes such as EMT. We show that this crosstalk originates downstream of the SMAD proteins phosphorylation. Using mathematical models, we predict, and experimentally verify, that promiscuous interactions between SMAD proteins provide the mechanism for the observed crosstalk. Our findings challenge the canonical models, suggesting an active role for mediator proteins in determining biological responses.

The adaptive immune system continuously encounters antigens from a wide range of sources, including pathogens, self-tissues, tumors, and environmental agents. While extensive studies have characterized how lymphocytes respond to antigen binding, most experimental frameworks consider the antigenic environment as static. In reality, antigen levels can fluctuate dramatically across a wide range of temporal and spatial scales. In this review, we examine how the dynamics of antigen presentation, ranging from molecular binding events to organism-level exposure, affect T cell activation and fate. We discuss the cellular and molecular mechanisms that allow T cells to detect and respond to changes in antigen concentration over timescales from seconds to days. These include kinetic proofreading of TCR signaling, frequency-dependent decoding in intracellular signaling networks, and population-level feedback circuits involving effector and regulatory T cells. Theoretical and experimental evidence suggests that T cells are tuned not only to antigen quantity but also to its rate of change, with implications for tolerance, immune activation, and memory formation. We highlight how manipulating the dynamics of antigen exposure, such as through controlled vaccine delivery, can modulate immune responses and suggest that incorporating temporal features into immunological models may improve our understanding of immune decision-making and inform therapeutic strategies.

Abstract Cells interpret complex extracellular environments through signaling pathways that compress diverse ligand inputs into a limited set of intracellular mediators. How this information is encoded in transcriptional responses, and whether it can be decoded to infer the original signaling environment, remains unclear. Here, we systematically map BMP/TGFβ transcriptional responses using single-cell RNA sequencing across 48 ligand-concentration conditions, generating a multi-ligand dose-response atlas. We find that, for each ligand, transcriptional responses collapse onto a one-dimensional, highly coordinated program in which concentration modulates response amplitude without altering gene identity. Across ligands, we find a small number of distinct, ligand-dependent transcriptional programs. We leverage this structured encoding to define a quantitative perception score that captures pathway activity at single-cell resolution. Finally, we train a machine learning model to decode extracellular ligand concentrations from single-cell transcriptomes. The model further generalizes to in vivo intestinal epithelium, where it reconstructs the spatial BMP gradient profile. Together, these results reveal a low-dimensional transcriptional code in BMP/TGFβ signaling that constrains cellular responses and enables quantitative inference of extracellular signaling environments from gene expression.

T helper 17 and Regulatory T cells (Th17 and Treg, respectively) are two well-described lymphocyte subsets with opposing actions. The divergent fates of Th17 and Treg cells are accounted for, at least in part, by molecular antagonism that occurs between their respective specific transcription factors, RORγt and Foxp3. An imbalance between Th17 and Treg cells may lead to tissue inflammation and is associated with certain types of autoimmunity. In order to understand the heterogeneity and dynamics of the differentiation process, we studied Th17/Treg cell differentiation of naïve cells in vitro, using RORγtGFPFoxp3RFP dual-reporter mouse. Flow cytometry revealed the consistent emergence of a population of double positive RORγt+Foxp3+ (DP) cells during the early stages of Th17 cell differentiation. These DP cells are closely related to RORγt+ single positive (SPR) cells in terms of global gene expression. Nevertheless, for some genes, DP cells share an expression pattern with Foxp3+ single positive (SPF) Treg cells, most importantly by reducing IL17 levels. Using time-lapse microscopy, we could delineate the expression dynamics of RORγt and Foxp3 at a clonal level. While the RORγt expression level elevates early during differentiation, Foxp3 rises later and is more stable upon environmental changes. These distinct expression profiles are independent of each other. During differentiation and proliferation, individual cells transit between SPR, DP, and SPF states. Nevertheless, the differentiation of sister cells within a clone progeny was highly correlated. We further demonstrated that sorted SPR and DP populations were not significantly affected by changes in their environment, suggesting that the correlated fate decision emerged at early time points before the first division. Overall, this study provides the first quantitative analysis of differentiation dynamics during the generation of DP RORγt+Foxp3+ cells. Characterizing these dynamics and the differentiation trajectory could provide a profound understanding and be used to better define the distinct fates of T cells, critical mediators of the immune response.

m6A is the most widespread mRNA modification and is primarily implicated in controlling mRNA stability. Fundamental questions pertaining to m6A are the extent to which it is dynamically modulated within cells and across stimuli, and the forces underlying such modulation. Prior work has focused on investigating active mechanisms governing m6A levels, such as recruitment of m6A writers or erasers leading to either global or site-specific modulation. Here, we propose that changes in m6A levels across subcellular compartments and biological trajectories may result from passive changes in gene-level mRNA metabolism. To predict the intricate interdependencies between m6A levels, mRNA localization, and mRNA decay, we establish a differential model m6ADyn encompassing mRNA transcription, methylation, export, and m6A-dependent and independent degradation. We validate the predictions of m6ADyn in the context of intracellular m6A dynamics, where m6ADyn predicts associations between relative mRNA localization and m6A levels, which we experimentally confirm. We further explore m6ADyn predictions pertaining to changes in m6A levels upon controlled perturbations of mRNA metabolism, which we also experimentally confirm. Finally, we demonstrate the relevance of m6ADyn in the context of cellular heat stress response, where genes subjected to altered mRNA product and export also display predictable changes in m6A levels, consistent with m6ADyn predictions. Our findings establish a framework for dissecting m6A dynamics and suggest the role of passive dynamics in shaping m6A levels in mammalian systems.

Signaling pathways enable cells to coordinate collective behaviors by exchanging specific information. Many pathways utilize multiple ligand variants to activate the same intracellular signaling cascade, raising the question of how cells discriminate between these seemingly redundant signals. It has been shown that individual cells can discriminate between signals based on their induced level of activity, temporal dynamics or combinatorial effect. Here, we demonstrate that ligand discrimination could also emerge at the population level. Using mathematical models of ligand-receptor interactions, we examine how response heterogeneity at the population level can encode ligand identity. We introduce a local scaling metric to quantify how variation in pathway components affects the cellular response. Our results reveal that for pathways with dimeric receptors, and more significantly for heterodimeric receptors, biochemical parameters of the ligands control the resulting heterogeneity in the response of a population of cells. Furthermore, we show that the population-level heterogeneity encodes the enzymatic activity of the resulting receptor complex. This suggests a functional advantage for utilizing heterodimeric receptor complexes in pathways acting across a population of cells, such as the type I interferon pathway, which shows several of the characteristics of our model. This contrasts to juxtacrine pathways, such as Notch, that do not act at the population level and use a single component receptor. Overall, our findings highlight a novel mechanism by which receptor architecture enables cells to encode ligand-specific information through population-level heterogeneity, offering insights into immune regulation, tissue development, and synthetic biology.

During an immune response, cells are simultaneously exposed to multiple cytokine signals that collectively determine their phenotype. Transforming growth factor β (TGFβ) is a pleiotropic cytokine acting as a key regulator of T-cell differentiation with activating and suppressive effects on their immune function. Here we identify and analyze the cellular responses of CD4+ T cells to TGFβ across signaling contexts by analyzing the responses of T cells to multiple cytokine mixtures in the presence or absence of TGFβ. We found that TGFβ had a profound dominant effect independent of the presence of other cytokines, modulating the expression of more than 4,000 genes. In the presence of TGFβ, cells exhibit lower expression of translation-related and apoptosis-related genes, accompanied by increased survival of activated T cells. Notably, cells cultured in the presence of TGFβ were smaller in size while preserving their proliferative ability. Accordingly, we identified a dense network of transcription factors that were modulated by TGFβ, suggesting a core gene set connecting TGFβ signaling to the regulation of T-cell size. We found N-Myc to be at the center of this network, and we directly show that TGFβ regulates its gene expression level, protein level, and nuclear localization. Our work provides a system to study cell size control and demonstrate the profound effect of TGFβ in the modulation and regulation of T-cell properties, expanding its role beyond guiding their phenotype. ### Competing Interest Statement The authors have declared no competing interest.

Numerous methods have recently emerged for ordering single cells along developmental trajectories. However, accurate depiction of developmental dynamics can only be achieved after rescaling the trajectory according to the relative time spent at each developmental point. We formulate a model which estimates local cell densities and fluxes, and incorporates cell division and apoptosis rates, to infer the real-time dimension of the developmental trajectory. We validate the model using mathematical simulations and apply it to experimental high dimensional cytometry data obtained from the mouse thymus to construct the true time profile of the thymocyte developmental process. Our method can easily be implemented in any of the existing tools for trajectory inference.

Flow cytometry enables quantitative measurements of fluorescence in single cells. The technique was widely used for immunology to identify populations with different surface protein markers. More recently, the usage of flow cytometry has been extended to additional readouts, including intracellular proteins and fluorescent protein transgenes, and is widely utilized to study development, systems biology, microbiology, and many other fields. A common file format (FCS format, defined by International Society for Advancement of Cytometry (ISAC)) has been universally adopted, facilitating data exchange between different machines. A diverse spectrum of software packages have been developed for analysis of flow cytometry data. However, those are either 1) costly proprietary softwares, 2) open source packages with prerequisite installation of R or Python and sometimes require users to have experience in coding or 3) online tools that are limiting for analysis of large data sets. Here we present EasyFlow, an open source flow cytometry analysis GUI based on Matlab or Python, that can be installed and run locally cross-platform-ly (Windows and MacOS), without prerequisite user having previous knowledge on coding. The python version (EasyFlowQ) is also developed on a popular plotting framework (Matplotlib) and modern user interface (UI) toolkit (Qt), allowing more advanced users to customize and keep contributing to the software, as well as its tutorials. Overall, EasyFlow serves as a simple-to-use tool for inexperienced users with little coding experience to use locally, as well as a platform for advanced users to further customize for their own needs.

Many biological circuits comprise sets of protein variants that interact with one another in a many-to-many, or promiscuous, fashion. These architectures can provide powerful computational capabilities that are espe-cially critical in multicellular organisms. Understanding the principles of biochemical computations in these circuits could allow more precise control of cellular behaviors. However, these systems are inherently difficult to analyze, due to their large number of interacting molecular components, partial redundancies, and cell context dependence. Here, we discuss recent experimental and theoretical advances that are beginning to reveal how promiscuous circuits compute, what roles those computations play in natural biological con-texts, and how promiscuous architectures can be applied for the design of synthetic multicellular behaviors.

Cell-cell communication systems typically comprise families of ligand and receptor variants that function together in combinations. Pathway activation depends on the complex way in which ligands are presented extracellularly and receptors are expressed by the signal-receiving cell. To understand the combinatorial logic of such a system, we systematically measured pairwise bone morphogenetic protein (BMP) ligand interactions in cells with varying receptor expression. Ligands could be classified into equivalence groups based on their profile of positive and negative synergies with other ligands. These groups varied with receptor expression, explaining how ligands can functionally replace each other in one context but not another. Context-dependent combinatorial interactions could be explained by a biochemical model based on the competitive formation of alternative signaling complexes with distinct activities. Together, these results provide insights into the roles of BMP combinations in developmental and therapeutic contexts and establish a framework for analyzing other combinatorial, context-dependent signaling systems.

In multicellular organisms, secreted ligands selectively activate, or “address,” specific target cell populations to control cell fate decision-making and other processes. Key cell-cell communication pathways use multiple promiscuously interacting ligands and receptors, provoking the question of how addressing specificity can emerge from molecular promiscuity. To investigate this issue, we developed a general mathematical modeling framework based on the bone morphogenetic protein (BMP) pathway architecture. We find that promiscuously interacting ligand-receptor systems allow a small number of ligands, acting in combinations, to address a larger number of individual cell types, defined by their receptor expression profiles. Promiscuous systems outperform seemingly more specific one-to-one signaling architectures in addressing capability. Combinatorial addressing extends to groups of cell types, is robust to receptor expression noise, grows more powerful with increases in the number of receptor variants, and is maximized by specific biochemical parameter relationships. Together, these results identify design principles governing cellular addressing by ligand combinations.

Numerous methods have recently emerged for ordering single cells along developmental trajectories. However, accurate depiction of developmental dynamics can only be achieved after rescaling the trajectory according to the relative time spent at each developmental point. We formulate a model which estimates local cell densities and fluxes, and incorporates cell division and apoptosis rates, to infer the real time dimension of the developmental trajectory. We validate the model using mathematical simulations, and apply it on experimental high dimensional cytometry data obtained from the mouse thymus to construct the true time-profile of the thymocyte developmental process. Our method can easily be implemented in any of the existing tools for trajectory inference.

Nascent messenger RNA is endowed with a poly(A) tail that is subject to gradual deadenylation and subsequent degradation in the cytoplasm. Deadenylation and degradation rates are typically correlated, rendering it difficult to dissect the determinants governing each of these processes and the mechanistic basis of their coupling. Here we developed an approach that allows systematic, robust and multiplexed quantification of poly(A) tails in Saccharomyces cerevisiae . Our results suggest that mRNA deadenylation and degradation rates are decoupled during meiosis, and that transcript length is a major determinant of deadenylation rates and a key contributor to reshaping of poly(A) tail lengths. Meiosis-specific decoupling also leads to unique positive associations between poly(A) tail length and gene expression. The decoupling is associated with a focal localization pattern of the RNA degradation factor Xrn1, and can be phenocopied by Xrn1 deletion under nonmeiotic conditions. Importantly, the association of transcript length with deadenylation rates is conserved across eukaryotes. Our study uncovers a factor that shapes deadenylation rate and reveals a unique context in which degradation is decoupled from deadenylation.

Adenosine-to-inosine editing is catalyzed by ADAR1 at thousands of sites transcriptome-wide. Despite intense interest in ADAR1 from physiological, bioengineering, and therapeutic perspectives, the rules of ADAR1 substrate selection are poorly understood. Here, we used large-scale systematic probing of ∼2,000 synthetic constructs to explore the structure and sequence context determining editability. We uncover two structural layers determining the formation and propagation of A-to-I editing, independent of sequence. First, editing is robustly induced at fixed intervals of 35 bp upstream and 30 bp downstream of structural disruptions. Second, editing is symmetrically introduced on opposite sites on a double-stranded structure. Our findings suggest a recursive model for RNA editing, whereby the structural alteration induced by the editing at one site iteratively gives rise to the formation of an additional editing site at a fixed periodicity, serving as a basis for the propagation of editing along and across both strands of double-stranded RNA structures.

The bone morphogenetic protein (BMP) signaling pathway comprises multiple ligands and receptors that interact promiscuously with one another and typically appear in combinations. This feature is often explained in terms of redundancy and regulatory flexibility, but it has remained unclear what signal-processing capabilities it provides. Here, we show that the BMP pathway processes multi-ligand inputs using a specific repertoire of computations, including ratiometric sensing, balance detection, and imbalance detection. These computations operate on the relative levels of different ligands and can arise directly from competitive receptor-ligand interactions. Furthermore, cells can select different computations to perform on the same ligand combination through expression of alternative sets of receptor variants. These results provide a direct signal-processing role for promiscuous receptor-ligand interactions and establish operational principles for quantitatively controlling cells with BMP ligands. Similar principles could apply to other promiscuous signaling pathways.

Animal cells use a conserved repertoire of intercellular signaling pathways to communicate with one another. These pathways are well-studied from a molecular point of view. However, we often lack an “operational” understanding that would allow us to use these pathways to rationally control cellular behaviors. This requires knowing what dynamic input features each pathway perceives and how it processes those inputs to control downstream processes. To address these questions, researchers have begun to reconstitute signaling pathways in living cells, analyzing their dynamic responses to stimuli, and developing new functional representations of their behavior. Here we review important insights obtained through these new approaches, and discuss challenges and opportunities in understanding signaling pathways from an operational point of view.

As they proliferate, living cells undergo transitions between specific molecularly and developmentally distinct states. Despite the functional centrality of these transitions in multicellular organisms, it has remained challenging to determine which transitions occur and at what rates without perturbations and cell engineering. Here, we introduce kin correlation analysis (KCA) and show that quantitative cell-state transition dynamics can be inferred, without direct observation, from the clustering of cell states on pedigrees (lineage trees). Combining KCA with pedigrees obtained from time-lapse imaging and endpoint single-molecule RNA-fluorescence in situ hybridization (RNA-FISH) measurements of gene expression, we determined the cell-state transition network of mouse embryonic stem (ES) cells. This analysis revealed that mouse ES cells exhibit stochastic and reversible transitions along a linear chain of states ranging from 2C-like to epiblast-like. Our approach is broadly applicable and may be applied to systems with irreversible transitions and non-stationary dynamics, such as in cancer and development.