We constructed a computational model to describe NF-κB activation events in response to IKK activation by TNFα (Hoffmann et al., 2002). This model comprised a singular NF-κB species, four IκB isoforms IκBα/β/ε, and IKK. Synthesis and degradation of the IκBs and cellular localization and interactions for all components are calculated using a system of ordinary differential equations.
More recent studies by our group expanded the model to include numerically defined IKK activity (Werner et al., 2005) and inducible syntheses of IκBε and IκBβ (Kearns et al., 2006). The most recent version includes nfkb2 p100/IκBδ as a fourth IκB protein that is degraded in response to IKK1 activity to release NF-κB/RelA dimers (Basak et al., 2007).
A web interface for the MATLAB versions of these models is available at the link below.
Below are links to our code, tutorials, and supplemental data. Each download is listed with its relevant publication (listed in reverse chronological order).
Humoral immunity depends on efficient activation of B-cells and their subsequent differentiation to antibody secreting cells (ASCs). The transcription factor NFkB cRel is critical for B-cell proliferation, but incorporating its known regulatory interactions into a mathematical model of the ASC differentiation circuit prevented ASC generation in simulations. Indeed, experimental ectopic cRel expression blocked ASC differentiation by inhibiting the transcription factor Blimp1, and in wild type cells cRel was dynamically repressed during ASC differentiation by Blimp1 binding the Rel locus. Including this bi-stable circuit of mutual cRel-Blimp1 antagonism into a multi-scale model revealed that dynamic repression of cRel controls the switch from B-cell proliferation to ASC generation phases and hence the respective cell population dynamics. Our studies provide a mechanistic explanation of how dysregulation of this bi-stable circuit may result in pathologic B-cell population phenotypes and thus offer new avenues for diagnostic stratification and treatment.
Rapid antibody production in response to invading pathogens requires the dramatic expansion of pathogen-derived antigenspecific B lymphocyte populations. Whether B cell population dynamics
are based on stochastic competition between competing cell fates, as in the development of competence by the bacterium Bacillus subtilis, or on deterministic cell fate decisions that execute a predictable program, as during the development of the worm Caenorhabditis elegans, remains unclear. Here, we developed long-term live-cell microscopy of B cell population expansion and multiscale mechanistic computational modeling to characterize the role of molecular noise in determining phenotype heterogeneity. We show that the cell lineage trees underlying B cell population dynamics are mediated by a largely predictable decision-making process where the heterogeneity of cell proliferation and death decisions at any given
timepoint largely derives from nongenetic heterogeneity in the founder cells. This means that contrary to previous models, only a minority of genetically identical founder cells contribute the majority to the population response. We computationally predict and experimentally confirm nongenetic molecular determinants that are predictive of founder cells’ proliferative capacity. While founder cell heterogeneity may arise from different exposure histories, we show that it may also be due to the gradual accumulation of small amounts of intrinsic noise during the lineage differentiation process of hematopoietic stem cells to mature B cells. Our finding of the largely deterministic nature of B lymphocyte responses may provide opportunities for diagnostic and therapeutic development.
Understanding the functions of multi-cellular organs in terms of the molecular networks within each cell is an important step in the quest to predict phenotype from genotype. B lymphocyte population dynamics, which are predictive of immune response and vaccine effectiveness, are determined by individual cells undergoing division or death seemingly stochastically. Based on tracking single-cell time-lapse trajectories of hundreds of B cells, single-cell transcriptome and immunofluorescence analyses, we constructed an agent-based multi-modular computational model to simulate lymphocyte population dynamics in terms of the molecular networks that control NF-κB signaling, the cell cycle and apoptosis. Combining modeling and experimentation, we found that NF-κB cRel enforces the execution of a cellular decision between mutually exclusive fates by promoting survival in growing cells. But as cRel-deficiency causes growing B-cells to die at similar rates to non-growing cells, our analysis reveals that the phenomenological decision model of wild type cells is rooted in a biased race of cell fates. We show that a multi-scale modeling approach allows for prediction of dynamic organ-level physiology in terms of intra-cellular molecular networks.
When messenger RNA splicing occurs cotranscriptionally, the potential for kinetic control based on transcription dynamics is widely recognized. Indeed, perturbation studies have reported that when transcription kinetics are perturbed genetically or pharmacologically splice patterns may change. However, whether kinetic control is contributing to the control of splicing within the normal range of physiological conditions remains unknown. We examined if the kinetic determinants for co-transcriptional splicing (CTS) might be re- flected in the structure and expression patterns of the genome and epigenome. To identify and then quantitatively relate multiple, simultaneous CTS determinants, we constructed a scalable mathematical model of the kinetic interplay of RNA synthesis and CTS and parameterized it with diverse next generation sequencing (NGS) data. We thus found a variety of CTS determinants encoded in vertebrate genomes and epigenomes, and that these combine variously for different groups of genes such as housekeeping versus regulated genes. Together, our findings indicate that the kinetic basis of splicing is functionally and physiologically relevant, and may meaningfully inform the analysis of genomic and epigenomic data to provide insights that are missed when relying on statistical approaches alone.
The immune response is a concerted dynamic multi-cellular process. Upon infection, the dynamics of lymphocyte populations are an aggregate of molecular processes that determine the activation, division, and longevity of individual cells. The timing of these single-cell processes is remarkably widely distributed with some cells undergoing their third division while others undergo their first. High cell-to-cell variability and technical noise pose challenges for interpreting popular dye-dilution experiments objectively. It remains an unresolved challenge to avoid under- or over-interpretation of such data when phenotyping gene-targeted mouse models or patient samples. Here we develop and characterize a computational methodology to parameterize a cell population model in the context of noisy dye-dilution data. To enable objective interpretation of model fits, our method estimates fit sensitivity and redundancy by stochastically sampling the solution landscape, calculating parameter sensitivities, and clustering to determine the maximum-likelihood solution ranges. Our methodology accounts for both technical and biological variability by using a cell fluorescence model as an adaptor during population model fitting, resulting in improved fit accuracy without the need for ad hoc objective functions. We have incorporated our methodology into an integrated phenotyping tool, FlowMax, and used it to analyze B cells from two NFkB knockout mice with distinct phenotypes; we not only confirm previously published findings at a fraction of the expended effort and cost, but reveal a novel phenotype of nfkb1/p105/50 in limiting the proliferative capacity of B cells following B-cell receptor stimulation. In addition to complementing experimental work, FlowMax is suitable for high throughput analysis of dye dilution studies within clinical and pharmacological screens with objective and quantitative conclusions.
The steady states of cells affect their response to perturbation. Indeed, diagnostic markers for predicting the response to therapeutic perturbation are often based on steady state measurements. In spite of this, no method exists to systematically characterize the relationship between steady state and response. Mathematical models are established tools for studying cellular responses, but characterizing their relationship to the steady state requires that it have a parametric, or analytical, expression. For some models, this expression can be derived by the King-Altman method. However, King-Altman requires that no substrate act as an enzyme, and is therefore not applicable to most models of signal transduction. For this reason we developed py-substitution, a simple but general method for deriving analytical expressions for the steady state of any mass action model. Where the King-Altman method is applicable, we show that py-substitution yields an equivalent expression, and at comparable efficiency. We use py-substitution to study the relationship between steady state and sensitivity to the anti-cancer drug candidate, dulanermin (recombinant human TRAIL). First, we use py-substitution to derive an analytical expression for the steady state of a published model of TRAIL-induced apoptosis. We show that the amount of TRAIL required for cell death is sensitive to the steady state concentrations of procaspase 8 and its negative regulator, Bar, but not the other procaspase molecules. This suggests that activation of caspase 8 is a critical point in the death decision process. Finally, we show that changes in the threshold at which TRAIL results in cell death is not always equivalent to changes in the time of death, as is commonly assumed. Our work demonstrates that an analytical expression is a powerful tool for identifying steady state determinants of the cellular response to perturbation.