Workflow for the construction, calibration, and analysis of Boolean models of vaccine-induced immune responses combining prior knowledge and context-specific experimental data.

Schematic overview of the methodology used to integrate prior biological knowledge and context-specific experimental data into executable Boolean networks. The workflow illustrates successive steps to build a Boolean network combining existing knowledge from different sources and experimental data to construct a naïve Boolean network and calibrate it respectively, and to analyze it (both statically and dynamically after simulation) to gain mechanistic insights and identify and test perturbation targets. Boolean networks were visualized using the Living Map software (BIOVIA, Dassault Systèmes). Proteins are represented as blue rectangles and phenotypes or cellular populations as orange hexagons. The heatmap illustrates node states across successive network updates, with green indicating active (on/1) and red inactive (off/O) states.

Broad description of the innate immune response to the MVA vaccine.

Top 25 Gene Ontology (GO) biological process terms identified by static GO enrichment analysis of the 200 genes included in the naïve Boolean network, compared to the human genome as reference. Dark grey bars correspond to processes intrinsic to the three constituent signaling pathways used to build the network (Cytosolic DNA-sensing, apoptosis, and NF-κB signaling KEGG pathways), whereas lighter grey bars indicate broader processes emerging from the integration of these pathways. The x-axis represents the negative logarithm of the associated p-values, all of which are below 0.05.

Simulated MVA-induced cellular abundance trajectories compared with binarized experimental data.

(a) Simulated binary immune cell population abundances, (b) Experimentally measured and binarized cellular population abundances. In panels (a) and (b), heatmaps indicate node states across successive network updates or experimental time points, with green representing active/high (1) and red inactive/low (0) states. (c) Comparison of simulated and experimental binary states. The heatmap represents correctly predicted states in blue and incorrectly predicted states in dark grey. Light grey columns correspond to intermediary simulation steps that do not match an experimental sampling time point.

Perturbation of Boolean network trajectories in simulations of the triple mutant MVA ΔC6L/ΔK7R/ΔA46R compared with parental MVA.

(a) Detailed heatmap showing downregulated (white), stable (pink), and upregulated (dark red) nodes for immune cellular populations and selected innate immune marker genes. Perturbations correspond to the simulated deletion of the C6L, K7R, and A46R genes in the MVA genome relative to the parental MVA condition. (b) System-wide heatmap displaying downregulated, stable, and upregulated nodes for the same perturbations scenario compared with parental MVA.

Perturbations of Boolean network trajectories in simulations of the MVA ΔN2L and MVA A21L-F55A/T116A/T117A mutants and the corresponding quadruple mutants MVA ΔN2L/ΔC6L/ΔK7R/ΔA46R and MVA A21L-F55A/T116A/T117A/ΔC6L/ΔK7R/ΔA46R compared with parental MVA.

Detailed heatmaps showing downregulated (white), stable (pink), and upregulated (dark red) nodes for immune cellular populations and selected innate immune marker genes. Perturbations correspond to: (a) deletion of the N2Lgene; (b) deletion of N2L, C6L, K7R, and A46R genes; (c) alanine substitutions F55A/T116A/T117A in the A21L gene; and (d) the same alanine substitutions in the A21L gene combined with deletion of C6L, K7R, and A46R genes. All conditions are shown relative to the parental MVA simulation.

Comparison of simulated and experimental binary immune cell population dynamics following MVA and YF-17D vaccination.

Simulated binary immune cell population abundances ((a) and (d)), experimentally measured binarized cellular abundances ((b) and (e)), and comparison between simulated and experimental binary states ((c) and (f)) are shown for MVA- ((a), (b), and (c)) and YF-17D-induced ((d), (e) and (f)) immune responses. Heatmaps display node states as On (green, 1) or Off (red, 0) across network updates or experimental time points. In the comparison panels, correctly predicted states are shown in blue and incorrectly predicted states in dark grey. Light grey columns correspond to intermediary simulation steps that do not match experimental time points.

Static analysis of the MVA- and YF-17D-induced immune response Boolean networks.

(a-b) Determinative Power (DP) vs. Vertex Betweenness (VB), reflecting node connectivity, for the MVA-induced (black, a) and YF-17D-induced (yellow, b) immune responses. Nodes with high DP and VB values are annotated, (c-d) Average effectiveness vs. number of outgoing edges, for the MVA-induced (black, c) and YF-17D-induced (yellow, d) immune responses. Nodes combining high average effectivity and large numbers of outgoing edges are annotated, Minimal values for annotation are represented as a red line. (e) Average effectiveness of regulatory functions in both networks as a function of the number of input variables (k). (f) Highly effective signaling paths (with relative effectiveness (→) > 0.6) leading to dynamically impactful nodes (in bold), defined as nodes with the highest connectivity (DP vs. VB), and average effectiveness. Blue boxes highlight host targets of immunomodulatory MVA genes previously described in the literature.

Simulated strategies for improving MVA viral vaccine vectors.

Dynamic simulation of immune cellular populations and stress-related MAPK marker genes following in silico deletion of the K7R (a) and F17R (b) genes from the MVA viral genome. Heatmaps represent downregulated (white), stable (pink), and upregulated (dark red) cell abundance or gene expression compared with the parental MVA.

Calibrated Boolean network of the MVA-induced immune response.

The network backbone is composed of the following three pathways from the KEGG database: cytosolic DNA-sensing, apoptosis, and NF-κB signaling. The naïve network was calibrated with experimental preclinical data from 46, using ZhegAlCal 74. The calibrated network is visualized using Dassault Systèmes Living Map software and is also available for simulation, perturbation, and analysis on the CellCollective platform (MVA 3 pathways).

Dynamic Gene Ontology (GO) enrichment analysis for the unperturbed MVA-induced innate immune response Boolean network trajectory.

The y-axis indicates the top 10 enriched biological process GO terms, while the x-axis displays the negative logarithm of the associated p-values (all < 0.05).

Binary cellular abundances and innate immune marker gene expressions in perturbed and unperturbed MVA-induced immune response Boolean network simulations.

(a) Parental MVA simulation. (b) Triple mutant MVA ΔC6L/ΔK7R/ΔA46R simulation.

Dynamic Gene Ontology enrichment analysis of the perturbed ΔC6L/ΔK7R/ΔA46R MVA-induced innate immune response Boolean network trajectory.

The y-axis indicates the top 10 enriched biological process GO terms, and the x-axis displays the negative logarithm of the associated p-values (all < 0.05).

Calibrated Boolean networks of the MVA- and YF-17D-induced immune responses.

The network backbone is composed of the following six pathways from the KEGG database: cytosolic DNA-sensing, apoptosis, Toll-like receptor signaling, RIG-I-like receptor signaling, NF-κB signaling and mTOR. The shared consensus naïve network was calibrated with experimental preclinical data from 46 and 71, and using ZhegAlCal74. The calibrated networks are visualized using Dassault Systèmes Living Map software and are also available for simulation, perturbation, and analysis on the CellCollective platform. (a) MVA-induced response Boolean network (MVA 6 pathways). (b) YF-17D-induced response Boolean network (YF17D).

Dynamic Gene Ontology enrichment analysis of the unperturbed MVA- and YF-17D-induced immune response Boolean network trajectories.

(a) MVA-induced response Boolean network. (b) YF-17D-induced response Boolean network. The y-axis indicates the top 10 enriched biological process GO terms and the x-axis displays the negative logarithm of the associated p-values (all < 0.05).

Binary cellular abundances and innate immune marker gene expression in perturbed and unperturbed MVA-induced immune response Boolean network simulations.

(a) Parental MVA, (b) MVA ΔK7R MVA, and (c) MVA ΔF17R simulations.