Build a Multi-Agent AI Workflow for Biological Network Modeling, Protein Interactions, Metabolism, and Cell Signaling Simulation
In this tutorial, we build a multi-agent workflow for biological systems modeling and explore how different computational components work together inside one unified systems biology pipeline. We generate synthetic biological data, analyze gene regulatory structure, predict protein-protein interactions, optimize metabolic pathway activity, and simulate a dynamic cell signaling cascade, all within a Colab environment that remains practical and reproducible. We also use an OpenAI model to act as a principal investigator, synthesizing the outputs of all specialized agents into a single expert-style biological interpretation that connects regulation, interaction networks, metabolism, and signaling into a broader scientific story. Copy CodeCopiedUse a different Browser import sys, subprocess, pkgutil def _install_if_missing(packages): missing = [] for p in packages: import_name = p[“import”] if pkgutil.find_loader(import_name) is None: missing.append(p[“pip”]) if missing: print(“Installing:”, “, “.join(missing)) subprocess.check_call([sys.executable, “-m”, “pip”, “install”, “-q”] + missing) _install_if_missing([ {“pip”: “openai”, “import”: “openai”}, {“pip”: “numpy”, “import”: “numpy”}, {“pip”: “pandas”, “import”: “pandas”}, {“pip”: “matplotlib”, “import”: “matplotlib”}, {“pip”: “networkx”, “import”: “networkx”}, {“pip”: “scikit-learn”, “import”: “sklearn”}, ]) import os import json import math import textwrap import random import getpass from dataclasses import dataclass from typing import Dict, List, Tuple, Any import numpy as np import pandas as pd import matplotlib.pyplot as plt import networkx as nx from sklearn.linear_model import LogisticRegression from sklearn.metrics import roc_auc_score, average_precision_score from sklearn.model_selection import train_test_split from sklearn.preprocessing import StandardScaler from openai import OpenAI np.random.seed(42) random.seed(42) OPENAI_API_KEY = None try: from google.colab import userdata OPENAI_API_KEY = userdata.get(“OPENAI_API_KEY”) if OPENAI_API_KEY: print(“Loaded OPENAI_API_KEY from Colab Secrets.”) except Exception: pass if not OPENAI_API_KEY: try: OPENAI_API_KEY = getpass.getpass(“Enter OPENAI_API_KEY (hidden input): “).strip() except Exception: OPENAI_API_KEY = input(“Enter OPENAI_API_KEY: “).strip() os.environ[“OPENAI_API_KEY”] = OPENAI_API_KEY client = OpenAI(api_key=OPENAI_API_KEY) OPENAI_MODEL = “gpt-4o-mini” We prepare the Colab environment and make sure all required libraries are available before the workflow begins. We import the scientific computing, machine learning, graph analysis, plotting, and OpenAI libraries that support the full biological systems pipeline from start to finish. We also securely load the OpenAI API key either from Colab Secrets or hidden input, initialize the client, and define the model so the notebook is ready for later LLM-based synthesis. Copy CodeCopiedUse a different Browser def sigmoid(x): return 1 / (1 + np.exp(-x)) def pretty(title: str, body: str, width: int = 100): print(“n” + “=” * width) print(title) print(“=” * width) print(body) def safe_float(x): try: return float(x) except Exception: return None def generate_gene_regulatory_network(n_genes: int = 14, edge_prob: float = 0.18): genes = [f”G{i+1}” for i in range(n_genes)] W = np.zeros((n_genes, n_genes)) for i in range(n_genes): for j in range(n_genes): if i != j and np.random.rand() < edge_prob: W[i, j] = np.random.uniform(-1.5, 1.5) return genes, W def simulate_gene_expression(W: np.ndarray, n_steps: int = 70, noise: float = 0.10): n = W.shape[0] X = np.zeros((n_steps, n)) X[0] = np.random.uniform(0.2, 0.8, size=n) for t in range(1, n_steps): signal = X[t-1] @ W X[t] = sigmoid(signal + np.random.normal(0, noise, size=n)) return X def generate_protein_features(n_proteins: int = 40, feature_dim: int = 10): proteins = [f”P{i+1}” for i in range(n_proteins)] features = np.random.normal(size=(n_proteins, feature_dim)) families = np.random.randint(0, 5, size=n_proteins) localization = np.random.randint(0, 4, size=n_proteins) return proteins, features, families, localization def generate_ppi_dataset(proteins, features, families, localization): rows = [] n = len(proteins) hidden_w = np.random.normal(size=features.shape[1]) for i in range(n): for j in range(i + 1, n): fi, fj = features[i], features[j] sim = np.dot(fi, fj) / (np.linalg.norm(fi) * np.linalg.norm(fj) + 1e-8) fam_same = 1 if families[i] == families[j] else 0 loc_same = 1 if localization[i] == localization[j] else 0 feat = np.concatenate([ np.abs(fi – fj), fi * fj, [sim, fam_same, loc_same] ]) score = 1.4 * sim + 1.0 * fam_same + 0.8 * loc_same + 0.15 * np.dot((fi + fj) / 2, hidden_w) prob = sigmoid(score) y = 1 if np.random.rand() < prob else 0 rows.append((proteins[i], proteins[j], feat, y)) return rows def generate_metabolic_network(): metabolites = [“Glucose”, “Pyruvate”, “AcetylCoA”, “ATP”, “Biomass”, “Lactate”, “Ethanol”] reactions = [ {“name”: “R1_Glucose_Uptake”, “yield_biomass”: 0.0, “yield_atp”: 0.3, “substrate_cost”: 1.0, “oxygen_need”: 0.0}, {“name”: “R2_Glycolysis”, “yield_biomass”: 0.2, “yield_atp”: 1.6, “substrate_cost”: 0.7, “oxygen_need”: 0.0}, {“name”: “R3_TCA”, “yield_biomass”: 1.0, “yield_atp”: 2.4, “substrate_cost”: 0.8, “oxygen_need”: 1.4}, {“name”: “R4_Fermentation”, “yield_biomass”: 0.1, “yield_atp”: 0.9, “substrate_cost”: 0.4, “oxygen_need”: 0.0}, {“name”: “R5_Ethanol_Path”, “yield_biomass”: 0.15,”yield_atp”: 0.8, “substrate_cost”: 0.5, “oxygen_need”: 0.0}, {“name”: “R6_Biomass_Assembly”,”yield_biomass”: 1.3, “yield_atp”: -0.9,”substrate_cost”: 0.6, “oxygen_need”: 0.2}, ] return metabolites, reactions def simulate_cell_signaling(T=200, dt=0.05, ligand_level=1.2): t = np.arange(0, T * dt, dt) ligand = np.ones_like(t) * ligand_level receptor = np.zeros_like(t) kinase = np.zeros_like(t) tf = np.zeros_like(t) phosphatase = np.zeros_like(t) receptor[0] = 0.05 kinase[0] = 0.02 tf[0] = 0.01 phosphatase[0] = 0.30 for i in range(1, len(t)): dR = 1.6 * ligand[i-1] * (1 – receptor[i-1]) – 0.9 * receptor[i-1] dK = 1.8 * receptor[i-1] * (1 – kinase[i-1]) – 1.1 * phosphatase[i-1] * kinase[i-1] dTF = 1.4 * kinase[i-1] * (1 – tf[i-1]) – 0.55 * tf[i-1] dP = 0.2 + 0.5 * tf[i-1] – 0.4 * phosphatase[i-1] receptor[i] = np.clip(receptor[i-1] + dt * dR, 0, 1) kinase[i] = np.clip(kinase[i-1] + dt * dK, 0, 1) tf[i] = np.clip(tf[i-1] + dt * dTF, 0, 1) phosphatase[i] = np.clip(phosphatase[i-1] + dt * dP, 0, 1.5) return pd.DataFrame({ “time”: t, “ligand”: ligand, “receptor_active”: receptor, “kinase_active”: kinase, “tf_active”: tf, “phosphatase”: phosphatase, }) We define the main helper utilities and all synthetic data generation functions that power the notebook’s biological tasks. We create functions for gene regulatory network construction, gene expression simulation, protein feature generation, protein interaction dataset creation, metabolic network setup, and cell signaling dynamics, which together provide four distinct biological views for analysis. This snippet forms the computational backbone of the tutorial by creating the structured inputs that each specialized agent will later process and interpret. Copy CodeCopiedUse a different Browser @dataclass class AgentResult: name: str summary: Dict[str, Any] class GeneRegulatoryNetworkAgent: def run(self, genes, W, X) -> AgentResult: corr = np.corrcoef(X.T) inferred_edges = [] true_edges = [] n = len(genes) for i in range(n): for j in range(n): if i == j: continue if abs(corr[i, j]) > 0.35: inferred_edges.append((genes[i], genes[j], float(corr[i, j]))) if abs(W[i, j]) > 1e-8: true_edges.append((genes[i], genes[j], float(W[i, j]))) centrality_graph = nx.DiGraph() for gi in genes: centrality_graph.add_node(gi) for i in range(n): for j in range(n): if abs(W[i, j]) >




